> #------------------------------------------------------------------- > # 8. GAMMA FAMILY AND CENSORED SURVIVAL DATA > #------------------------------------------------------------------- > > #------------------------------------------------------------------- > # INSULATE > #------------------------------------------------------------------- > m<-matrix(scan("/glim/datasets/insulate"),ncol=2,byrow=T) Read 152 items > kv<-m[,1] > time<-m[,2] > par(mfrow=c(2,2)) > plot(kv,log(time),main="Figure 8.1 Insulate") > summary(glm(time~1,family=Gamma(link=log))) Call: glm(formula = time ~ 1, family = Gamma(link = log)) Deviance Residuals: Min 1Q Median 3Q Max -3.620 -2.711 -2.121 -1.201 4.464 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 5.1429 0.1112 46.27 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for Gamma family taken to be 0.9390196) Null deviance: 371.64 on 75 degrees of freedom Residual deviance: 391.08 on 75 degrees of freedom AIC: 720.35 Number of Fisher Scoring iterations: 10 Warning message: Algorithm did not converge in: (if (is.empty.model(mt)) glm.fit.null else glm.fit)(x = X, y = Y, > glm<-glm(time~kv,family=Gamma(link=log)) > summary(glm) Call: glm(formula = time ~ kv, family = Gamma(link = log)) Deviance Residuals: Min 1Q Median 3Q Max -2.9057 -1.2852 -0.6397 0.2386 2.5891 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 21.41186 1.61239 13.28 <2e-16 *** kv -0.55537 0.04845 -11.46 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for Gamma family taken to be 1.734156) Null deviance: 371.64 on 75 degrees of freedom Residual deviance: 131.64 on 74 degrees of freedom AIC: 610.09 Number of Fisher Scoring iterations: 6 > abline(21.4118567 ,-0.5553653 ) > chisq<-(371.6379 -131.6428 )/1.734156 > c(chisq,1-pchisq(chisq,1)) [1] 138.3930 0.0000 > kv2<-kv^2 > summary(glm(time~kv+kv2,family=Gamma(link=log))) Call: glm(formula = time ~ kv + kv2, family = Gamma(link = log)) Deviance Residuals: Min 1Q Median 3Q Max -2.8593 -1.2722 -0.6321 0.2511 2.6274 Coefficients: Estimate Std. Error t value Pr(>|t|) (Intercept) 16.548968 14.058441 1.177 0.243 kv -0.250501 0.867247 -0.289 0.774 kv2 -0.004731 0.013294 -0.356 0.723 (Dispersion parameter for Gamma family taken to be 1.715972) Null deviance: 371.64 on 75 degrees of freedom Residual deviance: 131.43 on 73 degrees of freedom AIC: 611.93 Number of Fisher Scoring iterations: 6 > chisq<-(131.6428 -131.4256 )/1.715972 > c(chisq,1-pchisq(chisq,1)) [1] 0.1265755 0.7220095 > u<-1-exp(-time/fitted(glm)) > hist(u, main="Figure 8.2 Uniform?") > u<-sort(u) > uhat<-((1:76)-0.5)/76 > plot(uhat,u,xaxs="i",xlim=c(0,1),yaxs="i",ylim=c(0,1), main="Figure 8.3 P-P plot") > abline(0,1) > k<-sqrt(76)*(max(abs(u-uhat))+0.5/76) > k [1] 1.265117 > #------------------------------------------------------------------- > # USING THE INVERSE LIKELIHOOD > #------------------------------------------------------------------- > ystar<-1/time > summary(glm(ystar~1,family=poisson, weight=time)) Call: glm(formula = ystar ~ 1, family = poisson, weights = time) Deviance Residuals: Min 1Q Median 3Q Max -6.2317 0.8169 1.8586 2.5016 3.4640 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -4.5906 0.1147 -40.03 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 371.64 on 75 degrees of freedom Residual deviance: 371.64 on 75 degrees of freedom AIC: 808.8 Number of Fisher Scoring iterations: 6 > summary(glm(ystar~kv,family=poisson, weight=time)) Call: glm(formula = ystar ~ kv, family = poisson, weights = time) Deviance Residuals: Min 1Q Median 3Q Max -2.5889 -0.2386 0.6398 1.2851 2.9057 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -21.41087 1.26313 -16.95 <2e-16 *** kv 0.55534 0.03797 14.63 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 371.64 on 75 degrees of freedom Residual deviance: 131.64 on 74 degrees of freedom AIC: 570.81 Number of Fisher Scoring iterations: 7 > > #------------------------------------------------------------------- > # CENSORED EXPONENTIAL DATA > #------------------------------------------------------------------- > censor<-(time<10)*1.0 > censor.time<-censor*time+(1-censor)*10 > censor.time [1] 5.79 10.00 10.00 10.00 10.00 10.00 10.00 10.00 10.00 10.00 10.00 10.00 [13] 10.00 10.00 10.00 10.00 10.00 10.00 7.74 0.40 10.00 9.88 10.00 10.00 [25] 2.75 0.79 10.00 3.91 0.27 0.69 10.00 10.00 10.00 10.00 0.96 4.15 [37] 0.19 0.78 8.01 10.00 7.35 6.50 8.27 10.00 10.00 3.16 4.85 2.78 [49] 4.67 1.31 10.00 10.00 10.00 1.97 0.59 2.58 1.69 2.71 10.00 0.35 [61] 0.99 3.99 3.67 2.07 0.96 5.35 2.90 10.00 0.47 0.73 1.40 0.74 [73] 0.39 1.13 0.09 2.38 > plot(kv,log(time),main="Figure 8.4 Censored data",pch=".") > points(kv,log(censor.time)) > abline(log(10),0,lty=2) > abline(21.4118567 ,-0.5553653 ) > ystar<-censor/censor.time > summary(glm(ystar~1,family=poisson, weight=censor.time)) Call: glm(formula = ystar ~ 1, family = poisson, weights = censor.time) Deviance Residuals: Min 1Q Median 3Q Max -1.3789 -1.3789 0.3919 1.4339 2.7458 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -2.3531 0.1512 -15.56 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 156.36 on 75 degrees of freedom Residual deviance: 156.36 on 75 degrees of freedom AIC: 284.12 Number of Fisher Scoring iterations: 4 > summary(glm(ystar~kv,family=poisson, weight=censor.time)) Call: glm(formula = ystar ~ kv, family = poisson, weights = censor.time) Deviance Residuals: Min 1Q Median 3Q Max -2.49644 -0.63626 -0.05923 0.87314 2.57357 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -17.57026 2.43519 -7.215 5.39e-13 *** kv 0.45567 0.06996 6.514 7.33e-11 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 156.36 on 75 degrees of freedom Residual deviance: 103.47 on 74 degrees of freedom AIC: 233.23 Number of Fisher Scoring iterations: 5 > abline(17.5702038 , -0.4556724 ,lty=3) > > #------------------------------------------------------------------- > # LEUKEMIA > #------------------------------------------------------------------- > m<-matrix(scan("/glim/datasets/leukemx"),ncol=3,byrow=T) Read 126 items > sample<-factor(m[,1]) > weeks<-m[,2] > censor<-m[,3] > rate<-censor/weeks > summary(glm(rate~1,family=poisson,weight=weeks)) Call: glm(formula = rate ~ 1, family = poisson, weights = weeks) Deviance Residuals: Min 1Q Median 3Q Max -1.9702 -0.9532 0.3814 0.9308 1.9737 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -2.8922 0.1821 -15.89 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 54.503 on 41 degrees of freedom Residual deviance: 54.503 on 41 degrees of freedom AIC: 209.54 Number of Fisher Scoring iterations: 4 > summary(glm(rate~sample,family=poisson,weight=weeks)) Call: glm(formula = rate ~ sample, family = poisson, weights = weeks) Deviance Residuals: Min 1Q Median 3Q Max -1.32472 -0.75468 0.07899 0.87500 1.59679 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -3.6861 0.3332 -11.062 < 2e-16 *** sample2 1.5266 0.3983 3.833 0.000127 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 54.503 on 41 degrees of freedom Residual deviance: 38.017 on 40 degrees of freedom AIC: 195.05 Number of Fisher Scoring iterations: 5 > summary(glm(rate~1,family=binomial(link=cloglog),weight=weeks)) Call: glm(formula = rate ~ 1, family = binomial(link = cloglog), weights = weeks) Deviance Residuals: Min 1Q Median 3Q Max -1.9984 -0.9669 0.3943 0.9738 2.4051 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -2.8638 0.1826 -15.69 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 61.507 on 41 degrees of freedom Residual deviance: 61.507 on 41 degrees of freedom AIC: 233.84 Number of Fisher Scoring iterations: 5 > summary(glm(rate~sample,family=binomial(link=cloglog),weight=weeks)) Call: glm(formula = rate ~ sample, family = binomial(link = cloglog), weights = weeks) Deviance Residuals: Min 1Q Median 3Q Max -1.33318 -0.78489 0.08413 0.92541 2.07821 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -3.6734 0.3314 -11.086 < 2e-16 *** sample2 1.5745 0.3968 3.968 7.26e-05 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 61.507 on 41 degrees of freedom Residual deviance: 43.968 on 40 degrees of freedom AIC: 218.3 Number of Fisher Scoring iterations: 5 > > #------------------------------------------------------------------- > # AML > #------------------------------------------------------------------- > m<-matrix(scan("/glim/datasets/aml1"),ncol=13,byrow=T) Read 533 items > age<-m[,1] > surv <-m[,2] > censor <-m[,3] > rig <-m[,4] > rigg <-m[,5] > rigm <-m[,6] > riga <-m[,7] > aig <-m[,8] > aigg <-m[,9] > aigm <-m[,10] > aiga <-m[,11] > wbc <-m[,12] > lim <-m[,13] > rate<-censor/surv > summary(glm(rate~1,family=poisson,weight=surv)) Call: glm(formula = rate ~ 1, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -1.40862 -1.23825 -0.03934 1.18073 2.42371 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -6.6253 0.2182 -30.36 <2e-16 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 62.464 on 40 degrees of freedom AIC: 298.36 Number of Fisher Scoring iterations: 8 > summary(glm(rate~age+rig+rigg+rigm+riga+aig+aigg+aigm+aiga+wbc+lim,family=poisson,weight=surv)) Call: glm(formula = rate ~ age + rig + rigg + rigm + riga + aig + aigg + aigm + aiga + wbc + lim, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -1.7141 -0.9284 -0.4293 0.7597 2.3012 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -8.9067958 1.8545778 -4.803 1.57e-06 *** age -0.0006792 0.0167612 -0.041 0.9677 rig -1.2047454 1.1368475 -1.060 0.2893 rigg 1.4157193 1.2718617 1.113 0.2657 rigm 1.4221263 1.4110664 1.008 0.3135 aig 0.0064558 0.0081110 0.796 0.4261 aigg -0.0090970 0.0095732 -0.950 0.3420 aigm -0.0075873 0.0096159 -0.789 0.4301 wbc 0.3238362 0.2449703 1.322 0.1862 lim 0.0385363 0.0177197 2.175 0.0296 * --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 47.367 on 31 degrees of freedom AIC: 301.26 Number of Fisher Scoring iterations: 8 > summary(glm(rate~age+rig+aig+wbc+lim,family=poisson,weight=surv)) Call: glm(formula = rate ~ age + rig + aig + wbc + lim, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -1.7880 -0.9932 -0.4810 0.9959 2.3536 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -9.6652166 1.8357514 -5.265 1.40e-07 *** age 0.0004796 0.0153542 0.031 0.9751 rig 0.0558690 0.1370609 0.408 0.6836 aig -0.0011211 0.0007799 -1.438 0.1506 wbc 0.4202790 0.2314635 1.816 0.0694 . lim 0.0457534 0.0179807 2.545 0.0109 * --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 49.594 on 35 degrees of freedom AIC: 295.49 Number of Fisher Scoring iterations: 8 > summary(glm(rate~age+rig+wbc+lim,family=poisson,weight=surv)) Call: glm(formula = rate ~ age + rig + wbc + lim, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -1.7125 -0.9379 -0.3671 1.0084 2.3721 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -7.646880 1.103779 -6.928 4.27e-12 *** age 0.005247 0.015207 0.345 0.7301 rig -0.113440 0.071873 -1.578 0.1145 wbc 0.096911 0.052581 1.843 0.0653 . lim 0.025347 0.010719 2.365 0.0180 * --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 51.569 on 36 degrees of freedom AIC: 295.46 Number of Fisher Scoring iterations: 8 > summary(glm(rate~age+aig+wbc+lim,family=poisson,weight=surv)) Call: glm(formula = rate ~ age + aig + wbc + lim, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -1.6882 -1.0172 -0.4653 1.0144 2.3699 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -9.0517235 1.0155399 -8.913 < 2e-16 *** age 0.0015538 0.0150843 0.103 0.917954 aig -0.0008504 0.0004093 -2.078 0.037742 * wbc 0.3443430 0.1375585 2.503 0.012306 * lim 0.0402544 0.0117720 3.420 0.000627 *** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 49.756 on 36 degrees of freedom AIC: 293.65 Number of Fisher Scoring iterations: 8 > summary(glm(rate~age+wbc+lim,family=poisson,weight=surv)) Call: glm(formula = rate ~ age + wbc + lim, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -2.0493 -0.9902 -0.2009 0.9577 2.2861 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -8.724571 0.957154 -9.115 < 2e-16 *** age 0.008199 0.015837 0.518 0.60467 wbc 0.082416 0.053246 1.548 0.12167 lim 0.029815 0.010471 2.847 0.00441 ** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 54.525 on 37 degrees of freedom AIC: 296.42 Number of Fisher Scoring iterations: 8 > summary(glm(rate~wbc+lim,family=poisson,weight=surv)) Call: glm(formula = rate ~ wbc + lim, family = poisson, weights = surv) Deviance Residuals: Min 1Q Median 3Q Max -2.0178 -1.0181 -0.2374 0.9886 2.3067 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -8.41166 0.71287 -11.800 < 2e-16 *** wbc 0.08859 0.05090 1.740 0.08178 . lim 0.02951 0.01030 2.865 0.00417 ** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for poisson family taken to be 1) Null deviance: 62.464 on 40 degrees of freedom Residual deviance: 54.790 on 38 degrees of freedom AIC: 294.68 Number of Fisher Scoring iterations: 8 > > #------------------------------------------------------------------- > # MYELOMA > #------------------------------------------------------------------- > m<-matrix(scan("/glim/datasets/myeloma"),ncol=7,byrow=T) Read 455 items > censor<-m[,1] > time<-m[,2] > urea<-m[,3] > hemo<-m[,4] > age<-m[,5] > sex<-factor(m[,6]) > calcium<-m[,7] > rate<-censor/time > summary(glm(rate~urea+hemo+age+sex+calcium,family=binomial(link=cloglog),weight=time)) Call: glm(formula = rate ~ urea + hemo + age + sex + calcium, family = binomial(link = cloglog), weights = time) Deviance Residuals: Min 1Q Median 3Q Max -1.6634 -0.7624 0.2152 1.0959 2.0372 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -4.56543 1.63018 -2.801 0.00510 ** urea 1.70607 0.56521 3.018 0.00254 ** hemo -0.13032 0.06029 -2.161 0.03066 * age -0.02138 0.01553 -1.377 0.16863 sex1 -0.06088 0.31040 -0.196 0.84450 calcium 0.14304 0.09675 1.478 0.13927 --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 88.101 on 64 degrees of freedom Residual deviance: 72.953 on 59 degrees of freedom AIC: 425.62 Number of Fisher Scoring iterations: 6 > summary(glm(rate~urea,family=binomial(link=cloglog),weight=time)) Call: glm(formula = rate ~ urea, family = binomial(link = cloglog), weights = time) Deviance Residuals: Min 1Q Median 3Q Max -2.2813 -0.8227 0.3305 0.9907 2.2296 Coefficients: Estimate Std. Error z value Pr(>|z|) (Intercept) -5.7624 0.7794 -7.394 1.43e-13 *** urea 1.7114 0.5503 3.110 0.00187 ** --- Signif. codes: 0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 (Dispersion parameter for binomial family taken to be 1) Null deviance: 88.1 on 64 degrees of freedom Residual deviance: 79.8 on 63 degrees of freedom AIC: 424.47 Number of Fisher Scoring iterations: 5 > >