1根据以下表数据求生成函数在的Kaplan-Meier估计并用R

1.根据以下表数据求生成函数在t=204的Kaplan-Meier估计并用R绘出生存函数曲线。其中变量Xi=min（Ti，Ci）,1是第一次检查的时间，G是审查时间，示性函数可表示为：当第i个观察值缺失时为0;否则为1.Xi21418415070161412101323020484363869加11 11 1 1 0 1 0 1 1101分析：根据以上数据来求生成函数在t=204的Kaplan-Meier估计并用R绘出生存函数曲线。步骤：首先对观察死亡的时间进行排序：t⑴<...<3,di为在时刻如死亡的个体数；n为在时刻如面临危险的个体数，则在t⑴时刻危险函数的估计为：*udi/ni,从而得到生存函数在t时刻的Kaplan-Meier估计为：S(t)个」小式1一?)曰」st。一5)

t(i)— t(i)— n结论：经R语言计算得到：生成函数在 t=204的Kaplan-Meier估计为：0.1702381。R语言程序：cancer<-data.frame(times=c(214,184,150,70,16,141,210,132,30,204,84,36,38,69),status=c(1,1,1,1,1,1,0,1,0,1,1,1,0,1));with(cancer,Surv(times,status))summary(survfit(Surv(times,status),data=cancer))fit1<-survfit(Surv(times,status),data=cancer,conf.type="none")plot(fit1,main="KMforCancer",xlab="Time(days)",ylab="survivalrate")Call:survfit(formula=Surv(times,status),data=cancer)timen.riskn.eventsurvivalstd.errlower95%CIupper95%CITime(days)Time(days)161410.9290.06880.80301.000361210.8510.09730.68031.000691010.7660.11910.56481.00070910.6810.13290.46460.99884810.5960.14090.37480.947132710.5110.14420.29360.888141610.4260.14310.22020.823150510.3400.13750.15430.751184410.2550.12680.09650.676204310.1700.10940.04830.600214110.000NANANA(1-1/14)*(1-1/12)*(1-1/10)*(1-1/9)*(1-1/8)*(1-1/7)*(1-1/6)*(1-1/5)*(1-1/4)*(1-1/3)[1]0.1702381KMforCancerOrorlawvIwrus8060402.0OrorlawvIwrus8060402.0oHuoo5di2.根据一下数据求生存函数在t=1408的Kaplan-Meier估计。其中变量Xi=min(Ti,G),工是第一次检查的时间，G是审查时间，示性函数瓦表示为：当第i个观察值缺失时为0;否则为1.Xi346,141,296,1953,1375,822,2052,8361910,419,107,570,312加1,1,1, 0, 1, 1, 0, 1, 0, 1,1,1,1Xi1818,364,401,1645,330,1540,688,1309,505,1378,1446,86,1696i0,1,1, 0, 1, 0, 1, 0, 1,0, 0,1,0Xi2177,1968,1889,173,2070,1972,1897,2022,1879,1726,180Si0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0Xi615,1408,1763,1684,1576,1572,498,1585,1493,950,1242,119&i1, 1, 0, 0, 0, 0, 1, 0, 0,1, 0, 0分析：根据以上数据求生存函数在t=1408的Kaplan-Meier估计。步骤：首先对观察死亡的时间进行排序：t⑴<...<〃)，di为在时刻如死亡的个体数；n为在时刻％)面临危险的个体数，则在t⑴时刻危险函数的估计为：，?=di/ni,从而得到生存函数在t时刻的Kaplan-Meier估计为:S(t)=i「(…7)=Lt,t(i-di)

t(i)- t(i)-n结论：经R语言计算得到：生成函数在t=1408的Kaplan-Meier估计为：0.6084354。R语言程序：foof<-data.frame(times=c(346,141,296,1953,1375,822,2052,836,1910,419,107,570,312,1818,364,401,1645,330,1540,688,1309,505,1378,1446,86,1699,2177,1968,1889,173,2070,1972,1897,2022,1879,1726,1807,615,1408,1763,1684,1576,1572,498,1585,1493,950,1242,1190),status=c(1,1,1,0,1,1,0,1,0,1,1,1,1,0,1,1,0,1,0,1,0,1,0,0,1,0,0,0,0,1,0,0,0,0,0,0,0,1,1,0,0,0,0,1,0,0,1,0,0));with(cancer,Surv(times,status))summary(survfit(Surv(times,status),data=foof))Call:survfit(formula=Surv(times,status),data=foof)timen.riskn.eventsurvivalstd.errlower95%CIupper95%CI864910.9800.02020.9411.0001074810.9590.02830.9051.0001414710.9390.03420.8741.0001734610.9180.03910.8450.9982964510.8980.04320.8170.9873124410.8780.04680.7900.9743304310.8570.05000.7650.9613464210.8370.05280.7390.9473644110.8160.05530.7150.9324014010.7960.05760.6910.9174193910.7760.05960.6670.9024983810.7550.06140.6440.8865053710.7350.06310.6210.8695703610.7140.06450.5980.8536153510.6940.06580.5760.8366883410.6730.06700.5540.8188223310.6530.06800.5330.8018363210.6330.06890.5110.7839503110.6120.06960.4900.7651375 27 1 0.590 0.0706 0.466 0.7461408 25 1 0.566 0.0716 0.442 0.725(1-1/49)*(1-1/48)*(1-1/47)*(1-1/46)*(1-1/45)*(1-1/44)*(1-1/40)*(1-1/39)*(1-1/38)*(1-1/37)*(1-1/36)*(1-1/35)*(1-1/34)*(1-1/33)*(1-1/32)*(1-1/31)*(1-1/27)*(1-1/25)[1]0.60843541.给定危险率为九(t)=Kt,求相应的生存函数。解：由生存函数、危险率以及在[0,t]内的累积危险函数之间的关系：t L丁(y)dy=-lnS(t)=「：(t)从而得到：-lnS(t)=o,(y)dy=osJds=t所以相应的生存函数为：S(t)=e,''.给定危险率为九(t)=氐自口/。下+(t/。)，求在[0,t]内的累积危险函数。解：由生存函数、危险率以及在[0,t]内的累积危险函数之间的关系：JoMy)dy=-inS(t)=A(t)(累积危险函数)从而得到：A(t)=；(y)dy=°,y4/二'[1(y/o)']dy=ln[1(t/二)'].写出由Greenwood'Sormula,1926给出的在t时刻生存函数方差的估计式，同时根据以下表数据写出在时间t时生存函数的95%置信区间的R代码。其中Xi,2参考第一题的解释。Xi2828719586137747114022120176181155296i1 0 1 11 1 0 1 0 1 1101解:Greenwood'sformula,19第出的在t时刻生存函数方差的估计式为：V?(t)=S\t)% dit”ni(ni—di)其中：di为在时刻ti死亡的个体数；5为在时刻ti面临危险的个体数;S⑴=□t(i)<(1->?)=nt("(1-:)。在时间t时生存函数的95%置信区间的R代码为:cancer<-data.frame(times=c(28,287,195,86,137,74,71,140,22,120,176,181,155,29),status=c(1,0,1,1,1,1,0,1,0,1,1,1,0,1));fit4<-survfit(Surv(times,status),data=cancer,conf.type="plain")Call:survfit(formula=Surv(times,status),data=cancer,conf.type="plain")neventsmedian0.95LCL0.95UCL14 10 140 86 181.写出由Tsiatis,1981给出的在t时刻生存函数方差的估计式，同时根据以下表数据写出在时间t时生存函数的95%置信区间的R代码。其中Xi,。参考第一题的解释。Xi25296175541671502198611111281786i101110101111解:Tsiatis,1981给出的在t时刻生存函数方差的估计式为:Time(days)Time(days)A(t)=S2(t)\-2-

sni其中：di为在时刻ti死亡的个体数；n为在时刻ti面临危险的个体数;S(t)=rit<(「*)=口t/(1当。t(i)S t(i)& n在时间t时生存函数的95%置信区间的R代码为：cancer<-data.frame(times=c(252,96,175,54,167,150,219,86,1,111,128,178),status=c(1,0,1,1,1,0,1,0,1,1,1,1));fit6<-survfit(Surv(times,status),data=cancer)plot(fit6,main="KMcurveforOvarianCancerwith95%AalenCI",xlab="Time(days)",ylab="survivalrate")KMcurveforOvarianCancerwith95%AalenCIOrorlawvIwrus6040o2.0HuooOrorlawvIwrus6040o2.0Huoo10.在生存分析中，设di为在时刻t死亡(失效)的个体数，Q为在时刻ti面临危险的个体数，当djni很小时，寻求生成函数的Kaplan-Meier估计与Nelson-Aalen估计之间的关系。解：tA(t)=-lnS(t)^0(y)dy国t)=ilt。-？)="t.《/)t(i)S t(i)MR二—lnS(t)=—lnn (1—4=—£ln(1—5)%£5 (ln(1x)之x对充分小的x)t(i)Jn加上 ni“A所以，当di/ni彳艮小时,生成函数的Kaplan-Meier估计与Nelson-Aalen估计近似的相等。.根据以下表数据求生成函数在t=204的Nelson-Aalen估计。其中变量Xi=min(Ti,G),1是第一次检查的时间，G是审查时间，示性函数5表示为：当第i个观察值缺失时为0;否则为1.Xi21418415070161412101323020484363869露11 11 1 1 0 1 0 1 1101分析：根据以上数据求生成函数在t=204的Nelson-Aalen估计。步骤：首先对观察死亡的时间进行排序：t⑴<…<5,di为在时刻如死亡的个体数；n为在时刻t⑴面临危险的个体数，则在【0,端】累积危险函数的估计为：A⑴=11dL,从而得到生存函数在t时刻的Nelson-Aalen估计为：S(t)=e-At)。结论：经R语言计算得到：生成函数在t=204的Nelson-Aalen估计为:0.2052082。R语言程序：cancer<-data.frame(times=c(214,184,150,70,16,141,210,132,30,204,84,36,38,69),status=c(1,1,1,1,1,1,0,1,0,1,1,1,0,1));with(cancer,Surv(times,status))fit6<-coxph(Surv(times,status)~1,data=cancer)summary(survfit(fit6))Call:survfit.coxph.null(object=fit6)timen.riskn.eventsurvivalstd.errlower95%CIupper95%CI161410.93110.06650.80941.000361210.85660.09400.69081.000691010.77510.11510.57941.00070910.69360.12860.48220.99884810.61210.13690.39490.949132710.53060.14080.31540.893141610.44920.14080.24300.830150510.36770.13670.17750.762184410.28640.12830.11900.689204310.20520.11460.06870.613214110.07550.08650.00800.713x<-1/14+1/12+1/10+1/9+1/8+1/7+1/6+1/5+1/4+1/3exp(-x)[1]0.2052082.根据以下表数据求生成函数在t=i57的Nelson-Aalen估计。其中变量Xi=min（［，G）,［是第一次检查的时间，G是审查时间，示性函数百表示为：当第i个观察值缺失时为0;否则为1.Xi134,171,119,155,28,112,169,125,34,162,162,118,96,47,7：露1,1, 1, 1,1,1,0, 0,1,0,0,1,1,1,1Xi150,40,161,14,161,157,157,148,153,157,157,157,41d1,1,0,1,0, 0,1, 1,1, 0, 0, 0,1Xi151,46,151,71,151,150,150,148,12,18,70,470,1,0,1,0, 0,0,0,1,1,1,1分析：根据数据求生成函数在t=157的Nelson-Aalen估计。步骤：首先对观察死亡的时间进行排序：t（i）<...<t（k）,di为在时刻如死亡的个体数；门为在时刻t⑴面临危险的个体数，则在【0,t⑴】累积危险函数的估计为：A（t）=£2，从而得到生存函数在t时刻的5父nNelson-Aalen估计为：S（t）=e-At）。结论：经R语言计算得到：生成函数在 t=157的Nelson-Aalen估计为：0.7012939。R语言程序：k<-data.frame（times=c（134,171,119,155,28,112,169,125,34,162,162,118,96,47,73,150,40,161,14,161,157,157,148,153,157,157,157,41,151,46,151,71,151,150,150,148,12,18,70,47),

status=c(1,1,1,1,1,1,0,0,1,0,0,1,1,1,1,1,1,0,1,0,0,1,1,1,0,0,0,1,0,1,0,1,0,0,0,0,1,1,1,1));coxph(cancer,Surv(times,status))summary(survfit(Surv(times,status),data=k))Call:survfit(formula=Surv(times,status),data=k)timen.riskn.eventsurvivalstd.errlower95%CIupper95%CI124010.9750.02470.9281.000143910.9500.03450.8851.000183810.9250.04160.8471.000283710.9000.04740.8120.998343610.8750.05230.7780.984403510.8500.05650.7460.968413410.8250.06010.7150.952463310.8000.06320.6850.934473220.7500.06850.6270.897703010.7250.07060.5990.877712910.7000.07250.5710.857732810.6750.07410.5440.837962710.6500.07540.5180.8161122610.6250.07650.4920.7951182510.6000.07750.4660.7731192410.5750.07820.4410.7511342210.5490.07890.4140.7271482110.5230.07930.3880.7041501910.4950.07980.3610.6791531310.4570.08220.3210.6501551210.4190.08370.2830.6201571110.3810.08430.2470.588171110.000NANANAx<-1/40+1/39+1/38+1/37+1/36+1/35+1/34+1/33+1/32+1/30+1/29+1/28+1/27+1/26+1/25+1/24+1/22+1/21+1/19+1/13+1/12+1/11exp(-x)[1]0.7012939

.证明：具有固定带宽入和高斯核的Nadaraya-Watsonfc滑估计是可微的。对Epanechnikov核来说又怎么样了？证明：因为高斯核处处可微，所以对于给定的带宽九，K)(・,Xi)也是可/u微的，从而Nadaraya-Watson光滑估计是可微的。但是Epanechnikov核函数在t=±1是不可微的。.证明：对于局部线性回归，我们有'、Njx-5儿为人。证明：记b(xT)=x(1, w(X0)=diag(KMx—%),…,K,(Xn—X0)),li(x0)=(b(x0)T(BTw(%)B)"BTwd)}。DT=(x1-％，…,Xn-%).N, 、，， 一N, T,_T,、_、」_T,、一Jx0uy(Xi-X0)li(x0)=.：y(Xi-X0)li(x0)=b(x°)(Bw(Xo)B)Bw(xJx0u=b(x0)/、，) 、一x0=(1,Xo) =-X0<1 )=b(x=b(x0)/、，) 、一x0=(1,Xo) =-X0<1 )X0=0.证明任意阶的的局部多项式回归(包括局部常数) ，我们有：b0(x0)=1。其中："(％)=工;保-%)工函).证明：记b(x)T=(1,Xi,...,xd),w(Xo)=diag(K九(Xi—Xo),...,K*Xn—Xo)),li(x0)=(b(x0)T(BTw(x0)B)」BTw(x0))i,DT=((x-x°)j,...,(xN-x0)j)。由二项式定理：(a-b)n=C：anC：an4b...C：bn令：AT=(C0(—X0)j,C：(—X0)j\..,Cjj,0「,0)d-j1则有：D=BA。所以:

N j T_T _1_T _bj(xo)=1^Xi-X0)L(X0)=b(X0)(Bw(X0)B)-Bw(%)D=b(Xo)T(BTw(X0)B)」BTw(X0)Bb(X0)T=b(X0)TA//I、， 、,d\Cj=(1,X0,...,X0)0//I、， 、,d\Cj=(1,X0,...,X0)0=Cj(_X0)j，Cj(-Xo)'Xo ,CjX=(5-%),所以当j=0时有：b0（X0）=1.证明k阶的局部多项式回归，我们有：bj（%）=0,其中j=l,...,k;N jbj（X0）="a（Xj-x°）li（X0）.证明：由12题证明过程可以知道,4（％）=0,其中j=1,...,k。.在软件包”MASS里面，DDT包含了15个从农药DDT的不同实验室里面测得的指标，这些指标是利用多种农药残留检测的。用核方法找出DDT的密度曲线。解：library（MASS）plot(density(DDT))

density.default(x=DDT)2.53.0 3.5 4.04.55.0WWCLSneDdensity.default(x=DDT)2.53.0 3.5 4.04.55.0WWCLSneDN=15Bandwidth=0.1114.在软件包"boot”里面有salinity文件，它包含从1972到1977年的3月、4月和5月在PamlicoSound与北卡罗莱纳州两个地方测得的平均每两周组成的水含盐量和河流排量。用平均权重的Nadaraya-Watsonkernel及局部线性回归比较回归函数 E(平均每两周水含盐量|平均每两周河流排量)的估计。解：library(boot)plot(salinity$dis,salinity$sal,main="Nadaraya-Watsonkernel")lines(ksmooth(salinity$dis,salinity$sal),col=3,)library(KernSmooth)

g<-dpill(salinity$dis,salinity$sal)lines(locpoly(salinity$dis,salinity$sal,bandwidth=g),col=6)legend(29,14,legend=c("Nadaraya-Watson","Locallinear"),col=c(3,6),lty=1)Nadaraya-Watsonkernelsalinity$dis.普通线性回归模型y=P0+Pix+s其中日0凡分别是未知的载距与斜率参数，&是模型误差服从均值为0,方差为。2的正态分布。现给定数据（Xiy1x,y（2Xn）yn,.,..则（线，性）回归模型可表示为y+£凶=1,…,n其中鸟是独立同分布均值为0,方差为。2的正态分布，根据数据采用最小二乘法对未知的载距与斜率参数 P°,P1进行估计解：设bo心分别为Po,Pi的最小二乘估计，则拟合值为y=bo+b〔Xi,i=1,…，n,令y=bo+"(为—x),b0=b0+b1x,则最小二乘的目标n n a函数为：Zri2=£(yi-yi)2达到最小i1 i4从而有：口n _—vZ[y-bo-bi(x-x)]2=o-：boi4n

* 2[yi-bo-bi(xi-x)]-0>:bii4运用导数知识可以得到:n n*nbo+biZ(x-x)=£yi昌 i刍*n _ n _ n _bo*£(x—x)+b£(x-x)2=Syi(K-x)i1 i1 i1n _因为Z(x-x)=o,i1n

工yi(xi-x)…b 角牛伊： £(x-x)2i苴* _bo=y「n工yG-x)〜r 所以：£(x-x)2i苴一一bo=y-bix17.设f(x)具有p阶有界导数且p阶导数在f(x)支撑的内点x是连续的，假设函数K(x)满足UupK(u)du<z,则当hT。时，我们有:f(xhu)K(u)du=J-=Op'、Ui(K)f⑴(x)hii❷i!o(hp)其中：Uj(K)=uiK(u)du"*,"^O证明:由泰勒公式可以得到：p」

(i)ii

'f(x)huf(xhu): i!.f(P)(x(hu))hp.p!这里<hu)介于0与hu之间，所以:这里，(h)=xUi(K)f(i)(x)hif(xhu)K(u)du=i'二 i!hpJf(p)(x(hu))-f(p)(x)]upK(u)du.：(h)注意到f⑼⑷是有界的，由控制收敛定理可以得到A(h)=0(hp),证毕18.利用由伍斯特心脏病发作研究(WHAS)分组队列1(1975年至1978年)的数据。采用Log-rank 检验零假设：男性和女性具有相同的生存曲线。(数据在附件里)分析：根据数据采用Log-rank检验：零假设：男性和女性具有相同的生存曲线；备择假设：男性和女性具有不相同的生存曲线；步骤:令t⑴是在整个组中第i个人死亡的时间，dli表示男性在t(i)时刻死亡的人数，d2i表示女性在3时刻死亡的人数，小为在时刻t⑴面临危险的男性的人数，n2i为在时刻t⑴面临危险的女性的人数，记:di=diid2i n=n.n2iniidii RinzdiSi-di)eli二 v1i= 2~ ~ni 4(4-1)Oi-Ei八idii-^i Vi-¥则Log-rank 检验统计量为：经R语言计算得到：p=6.i3e-05结论：p=6.i3e-05 远远小于0.05,这意味着我们不应该接受零假设，也就是说，男性和女性具有不相同的生存曲线。R语言程序：首先把Excel数据以CSV(逗号分隔)(*csv)保存到F:盘。Library(splines)a<-read.csv("F:〃whas.csv",h=T)survdiff(Surv(cpk,chf)~sex,data=a)Call:survdiff(formula=Surv(cpk,chf)~sex,data=a)NObservedExpected(O-E)A2/E(O-E)A2/Vsex=0287 i02i28.4 5.43i6.isex=ii94 94 67.6 i0.30i6.iChisq=i6.ionidegreesoffreedom,p=6.i3e-0519.根据以下表数据求生成函数在t=170的Kaplan-Meier估计并用R绘出生存函数曲线。其中变量Xi=min（［，G）,［是第一次检查的时间，G是审查时间，示性函数可表示为：当第i个观察值缺失时为0;否则为1.Xi187,152,46,103,37,170,64,182,180,176,175,64,42,175,14!加0,1,1, 1,1,1,1,0,0, 0,0,1,1,0,1Xi158,169,169,22,168,157,89,165,14,161,161,1596i1,0, 0,0,0,1, 1,1,1,0,0,0分析：根据以上数据求生成函数在t=170的Kaplan-Meier估计并用R绘出生存函数曲线。步骤：首先对观察死亡的时间进行排序：t（1）<...<t（k）,di为在时刻如死亡的个体数；n为在时刻如面临危险的个体数，则在t⑴时刻危险函数的估计为：/udi/ni,从而得到生存函数在t时刻的Kaplan-Meier估计为：S（t）t.：t（1-？）t.：t（1-与，（i）二 yd二 「结论：经R语言计算得到：生成函数在 t=170的Kaplan-Meier估计为:0.4412121。R语言程序：cancer<-data.frame(times=c(187,152,46,103,37,170,64,182,180,176,175,64,42,175,149,158,169,169,22,168,157,89,165,14,161,161,159),status=c(0,1,1,1,1,1,1,0,0,0,0,1,1,0,1,1,0,0,0,0,1,1,1,1,0,0,0));with(cancer,Surv(times,status))summary(survfit(Surv(times,status),data=cancer))fit1<-survfit(Surv(times,status),data=cancer,conf.type="none")plot(fit1,main="KMforCancer",xlab="Time(days)",ylab="survivalrate")Call:survfit(formula=Surv(times,status),data=cancer)timen.riskn.eventsurvivalstd.errlower95%CIupper95%CI142710.9630.03630.8941.000372510.9240.05140.8291.000422410.8860.06200.7721.000462310.8470.07030.7200.997642220.7700.08230.6250.950892010.7320.08680.5800.9231031910.6930.09030.5370.8951491810.6550.09320.4950.8651521710.6160.09530.4550.8351571610.5780.09680.4160.8021581510.5390.09770.3780.7691651110.4900.10040.3280.732170710.4200.10780.2540