kaplan-meier.txt

kaplan-Meier估计,在可靠性数学中有一定的帮助

#Kaplan-Meier估计: 
KM<-function(x,v)
{n=length(x);if(n!=length(v))return("1")   #定义数据x的长度与指示值v的长度一样
 k=n-sum(v)        #系统总数减去截尾值的个数
 su=c(n:1)         #在失效时刻瞬间存活的系统数
 tf=c(1:k-1:k)     #失效时刻
 m=c(1:k-0:(k-1))  #失效数
 s=tf              #系统总数减去失效数和截尾系统数之和
 y=order(x)        #把数据x排序
 j=0             
 for(i in 1:n)
 {if(v[y[i]]==0){j=j+1;s[j]=su[i];tf[j]=x[i]}
}
w=c(1:k-1:k);
for(i in 1:(k-1))
{if(tf[i]==tf[i+1])w[i]=1}
l=sum(w)
ft=c(1:(k-l)-1:(k-l))  #存放失效时刻
so=ft
mo=ft
j=0
for(i in 1:k)
{if(w[i]==0)
   {j=j+1;so[j]=s[i];ft[j]=tf[i];mo[j]=m[i]}
 else m[i+1]=m[i]+m[i+1]}
Rel=c(1:(k-l)-1:(k-l))
VarofR=c(1:(k-l)-1:(k-l))
Rel[1]=1-mo[1]/so[1]                   #可靠性估计(KM估计)
VarofR=mo[1]^2/(so[1]*(so[1]-mo[1]))   #KM估计的渐进方差
for(i in 2:(k-l))
{
 Rel[i]=Rel[i-1]*(1-(mo[i]/so[i]))
 VarofR[i]=VarofR[i-1]+(mo[i]^2/(so[i]*(so[i]-mo[i])))
}
for(i in 1:(k-l))
VarofR[i]=Rel[i]*sqrt(VarofR[i]) #可靠性估计的标准差函数
return(ft,Rel,VarofR,so,mo)}
x=c(69,176,196,208,215,233,289,300,384,390,393,401,452,567,617,718,783,783,806,1000,1022)
v=c(0,0,1,0,0,0,0,0,0,0,1,0,0,0,1,1,0,0,0,1,0)
KM(x,v)

OUTPUT:
ft
 [1]   69  176  208  215  233  289  300  384  390  401  452  567  783  806 1022

Rel
 [1] 0.9523810 0.9047619 0.8544974 0.8042328 0.7539683 0.7037037 0.6534392
 [8] 0.6031746 0.5529101 0.4976190 0.4423280 0.3870370 0.1935185 0.1290123
[15] 0.0000000

VarofR
 [1] 0.04647143 0.06405645 0.07775694 0.08794121 0.09573806 0.10169808
 [7] 0.10613118 0.10922343 0.11108686 0.11290271 0.11308907 0.11165412
[13] 0.14778826 0.11171977        NaN

so
 [1] 21 20 18 17 16 15 14 13 12 10  9  8  4  3  1

mo
 [1] 1 1 1 1 1 1 1 1 1 1 1 1 2 1 1