Survival analysis of water quality data
Prepared by Byron and Anshul
Version – 16 May 2019
1 Dependent Variable: Functionality
This section contains descriptive statistics related to the TTC measure in the data.
TTC is level of fecal contamination. More TTC means more contamination in the sampled water.
Here are summary statistics for TTC:
## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 0.00 0.00 10.25 41.29 107.25 125.00 20## db2$TTC
## n missing distinct Info Mean Gmd .05 .10
## 330 0 78 0.953 41.34 53.07 0.00 0.00
## .25 .50 .75 .90 .95
## 0.00 10.25 107.25 125.00 125.00
##
## lowest : 0.0 0.5 1.0 1.5 2.0, highest: 99.0 110.0 113.0 121.0 125.0Below is a histogram of the TTC variable, which ranges from 0 to 125, with 0 meaning no contamination and 125 meaning 125 or higher level of contamination.
As you can see, the vast majority of the readings are either 0 or 125 and over.
If we break down the TTC variable into categories of any contamination (TTC > 0) or no contamination (TTC = 0), we have the following distribution:
| Any Contamination? | Frequency | Total | |
|---|---|---|---|
| Count | NA | ||
| 0 |
107 30.2 % |
0 0 % |
107 30.2 % |
| 1 |
227 64.1 % |
0 0 % |
227 64.1 % |
| NA |
20 5.6 % |
0 0 % |
20 5.6 % |
| Total |
354 100 % |
0 100 % |
354 100 % |
107 water sources have no contamination. 227 have contamination that is >0.
2 Independent Variables
2.1 Water Source Management Type
##
## Count <NA>
## Borehole handpump 144 0
## Borehole mech pump 22 0
## Dug well open 46 0
## Dug well handpump 71 0
## Spring protected 40 0
## Spring unprotected 10 0
## Rain harvest tank 19 0
## <NA> 2 02.2 Topography
##
## Count <NA>
## Low Lying 153 0
## Elevated 12 0
## Flat 187 0
## <NA> 2 03 Survival analysis
3.1 Preparation and analysis
Prepare the data:
# Prepare the data
db2$event <- NA
db2$event[db2$TTC<125] <- 1
db2$event[db2$TTC>=125] <- 0
# with(db, table(TTC,event, useNA='always'))
Sobj <- with(db2, Surv(TTC, event))Summarize survival data:
Sobj_fit <- survfit(Sobj ~ 1, data=db2)
# head(Sobj, 50)
# head(db2[c("TTC")], n=50)
summary(Sobj_fit, times = (0:130))## Call: survfit(formula = Sobj ~ 1, data = db2)
##
## time n.risk n.event survival std.err lower 95% CI upper 95% CI
## 0 330 106 0.679 0.0257 0.630 0.731
## 1 215 12 0.642 0.0264 0.593 0.696
## 2 208 9 0.615 0.0268 0.565 0.670
## 3 201 5 0.600 0.0270 0.549 0.655
## 4 197 2 0.594 0.0270 0.543 0.649
## 5 193 3 0.585 0.0271 0.534 0.641
## 6 190 5 0.570 0.0273 0.519 0.626
## 7 185 7 0.548 0.0274 0.497 0.605
## 8 175 8 0.524 0.0275 0.473 0.581
## 9 173 2 0.518 0.0275 0.467 0.575
## 10 166 6 0.500 0.0275 0.449 0.557
## 11 163 4 0.488 0.0275 0.437 0.545
## 12 159 2 0.482 0.0275 0.431 0.539
## 13 158 3 0.473 0.0275 0.422 0.530
## 14 156 1 0.470 0.0275 0.419 0.527
## 15 154 4 0.458 0.0274 0.407 0.515
## 16 148 3 0.448 0.0274 0.398 0.505
## 17 147 3 0.439 0.0273 0.389 0.496
## 18 144 2 0.433 0.0273 0.383 0.490
## 19 143 0 0.433 0.0273 0.383 0.490
## 20 142 3 0.424 0.0272 0.374 0.481
## 21 140 2 0.418 0.0272 0.368 0.475
## 22 137 1 0.415 0.0271 0.365 0.472
## 23 136 1 0.412 0.0271 0.362 0.469
## 24 136 2 0.406 0.0270 0.356 0.463
## 25 133 4 0.394 0.0269 0.345 0.450
## 26 126 5 0.379 0.0267 0.330 0.435
## 27 125 0 0.379 0.0267 0.330 0.435
## 28 125 0 0.379 0.0267 0.330 0.435
## 29 125 1 0.376 0.0267 0.327 0.432
## 30 123 3 0.367 0.0265 0.318 0.423
## 31 121 0 0.367 0.0265 0.318 0.423
## 32 121 0 0.367 0.0265 0.318 0.423
## 33 121 0 0.367 0.0265 0.318 0.423
## 34 121 0 0.367 0.0265 0.318 0.423
## 35 119 2 0.361 0.0264 0.312 0.416
## 36 119 0 0.361 0.0264 0.312 0.416
## 37 119 1 0.358 0.0264 0.309 0.413
## 38 117 2 0.352 0.0263 0.304 0.407
## 39 116 0 0.352 0.0263 0.304 0.407
## 40 116 1 0.348 0.0262 0.301 0.404
## 41 113 2 0.342 0.0261 0.295 0.398
## 42 113 2 0.336 0.0260 0.289 0.391
## 43 111 0 0.336 0.0260 0.289 0.391
## 44 111 2 0.330 0.0259 0.283 0.385
## 45 109 1 0.327 0.0258 0.280 0.382
## 46 108 1 0.324 0.0258 0.277 0.379
## 47 107 1 0.321 0.0257 0.275 0.376
## 48 106 0 0.321 0.0257 0.275 0.376
## 49 106 0 0.321 0.0257 0.275 0.376
## 50 106 3 0.312 0.0255 0.266 0.366
## 51 103 2 0.306 0.0254 0.260 0.360
## 52 101 0 0.306 0.0254 0.260 0.360
## 53 101 1 0.303 0.0253 0.257 0.357
## 54 100 2 0.297 0.0252 0.252 0.351
## 55 98 0 0.297 0.0252 0.252 0.351
## 56 98 2 0.291 0.0250 0.246 0.344
## 57 96 0 0.291 0.0250 0.246 0.344
## 58 95 1 0.288 0.0249 0.243 0.341
## 59 95 0 0.288 0.0249 0.243 0.341
## 60 95 1 0.285 0.0248 0.240 0.338
## 61 94 0 0.285 0.0248 0.240 0.338
## 62 94 0 0.285 0.0248 0.240 0.338
## 63 94 0 0.285 0.0248 0.240 0.338
## 64 94 0 0.285 0.0248 0.240 0.338
## 65 94 0 0.285 0.0248 0.240 0.338
## 66 94 0 0.285 0.0248 0.240 0.338
## 67 94 0 0.285 0.0248 0.240 0.338
## 68 94 1 0.282 0.0248 0.237 0.335
## 69 93 1 0.279 0.0247 0.234 0.332
## 70 92 0 0.279 0.0247 0.234 0.332
## 71 92 0 0.279 0.0247 0.234 0.332
## 72 92 0 0.279 0.0247 0.234 0.332
## 73 92 0 0.279 0.0247 0.234 0.332
## 74 92 0 0.279 0.0247 0.234 0.332
## 75 91 1 0.276 0.0246 0.232 0.328
## 76 91 0 0.276 0.0246 0.232 0.328
## 77 91 0 0.276 0.0246 0.232 0.328
## 78 91 0 0.276 0.0246 0.232 0.328
## 79 91 0 0.276 0.0246 0.232 0.328
## 80 91 0 0.276 0.0246 0.232 0.328
## 81 91 1 0.273 0.0245 0.229 0.325
## 82 90 0 0.273 0.0245 0.229 0.325
## 83 90 0 0.273 0.0245 0.229 0.325
## 84 90 0 0.273 0.0245 0.229 0.325
## 85 89 1 0.270 0.0244 0.226 0.322
## 86 89 0 0.270 0.0244 0.226 0.322
## 87 89 1 0.267 0.0243 0.223 0.319
## 88 87 1 0.264 0.0243 0.220 0.316
## 89 87 0 0.264 0.0243 0.220 0.316
## 90 86 1 0.261 0.0242 0.217 0.313
## 91 86 1 0.258 0.0241 0.214 0.309
## 92 84 1 0.255 0.0240 0.212 0.306
## 93 84 0 0.255 0.0240 0.212 0.306
## 94 84 0 0.255 0.0240 0.212 0.306
## 95 84 0 0.255 0.0240 0.212 0.306
## 96 84 0 0.255 0.0240 0.212 0.306
## 97 84 0 0.255 0.0240 0.212 0.306
## 98 84 0 0.255 0.0240 0.212 0.306
## 99 84 1 0.252 0.0239 0.209 0.303
## 100 83 0 0.252 0.0239 0.209 0.303
## 101 83 0 0.252 0.0239 0.209 0.303
## 102 83 0 0.252 0.0239 0.209 0.303
## 103 83 0 0.252 0.0239 0.209 0.303
## 104 83 0 0.252 0.0239 0.209 0.303
## 105 83 0 0.252 0.0239 0.209 0.303
## 106 83 0 0.252 0.0239 0.209 0.303
## 107 83 0 0.252 0.0239 0.209 0.303
## 108 83 0 0.252 0.0239 0.209 0.303
## 109 83 0 0.252 0.0239 0.209 0.303
## 110 83 1 0.248 0.0238 0.206 0.300
## 111 82 0 0.248 0.0238 0.206 0.300
## 112 82 0 0.248 0.0238 0.206 0.300
## 113 82 1 0.245 0.0237 0.203 0.297
## 114 81 0 0.245 0.0237 0.203 0.297
## 115 81 0 0.245 0.0237 0.203 0.297
## 116 81 0 0.245 0.0237 0.203 0.297
## 117 81 0 0.245 0.0237 0.203 0.297
## 118 81 0 0.245 0.0237 0.203 0.297
## 119 81 0 0.245 0.0237 0.203 0.297
## 120 81 0 0.245 0.0237 0.203 0.297
## 121 81 1 0.242 0.0236 0.200 0.293
## 122 80 0 0.242 0.0236 0.200 0.293
## 123 80 0 0.242 0.0236 0.200 0.293
## 124 80 0 0.242 0.0236 0.200 0.293
## 125 80 0 0.242 0.0236 0.200 0.293plot(Sobj_fit, xlab = "TTC", ylab = "Proportion continuing", main = "Baseline hazard curve")
Stratified by water source type:
Sobj_fit_Type <- survfit(Sobj ~ Type, data=db2)
plot(Sobj_fit_Type)
3.1.1 Cox proportional hazard model, few covariates
Sobj_cox <- coxph(Sobj ~ Type + Topo, data=db2)
summary(Sobj_cox)## Call:
## coxph(formula = Sobj ~ Type + Topo, data = db2)
##
## n= 330, number of events= 250
##
## coef exp(coef) se(coef) z Pr(>|z|)
## TypeBorehole mech pump -0.2396 0.7870 0.2685 -0.892 0.372221
## TypeDug well open -0.8797 0.4149 0.2504 -3.514 0.000442 ***
## TypeDug well handpump -0.2793 0.7563 0.1712 -1.632 0.102773
## TypeSpring protected -0.1872 0.8293 0.2110 -0.887 0.374892
## TypeSpring unprotected -1.7782 0.1689 0.7212 -2.466 0.013681 *
## TypeRain harvest tank 0.3808 1.4634 0.2553 1.492 0.135819
## TopoElevated 0.9329 2.5418 0.3202 2.913 0.003577 **
## TopoFlat 0.6520 1.9193 0.1434 4.548 5.42e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## TypeBorehole mech pump 0.7870 1.2707 0.4650 1.3320
## TypeDug well open 0.4149 2.4102 0.2540 0.6777
## TypeDug well handpump 0.7563 1.3222 0.5408 1.0578
## TypeSpring protected 0.8293 1.2059 0.5484 1.2539
## TypeSpring unprotected 0.1689 5.9190 0.0411 0.6945
## TypeRain harvest tank 1.4634 0.6833 0.8873 2.4135
## TopoElevated 2.5418 0.3934 1.3570 4.7612
## TopoFlat 1.9193 0.5210 1.4492 2.5420
##
## Concordance= 0.666 (se = 0.019 )
## Rsquare= 0.201 (max possible= 1 )
## Likelihood ratio test= 74.21 on 8 df, p=7e-13
## Wald test = 61.27 on 8 df, p=3e-10
## Score (logrank) test = 69.43 on 8 df, p=6e-12tab_model(Sobj_cox)| Sobj | |||
|---|---|---|---|
| Predictors | Estimates | CI | p |
| Borehole mech pump | 0.79 | 0.46 – 1.33 | 0.372 |
| Dug well open | 0.41 | 0.25 – 0.68 | <0.001 |
| Dug well handpump | 0.76 | 0.54 – 1.06 | 0.103 |
| Spring protected | 0.83 | 0.55 – 1.25 | 0.375 |
| Spring unprotected | 0.17 | 0.04 – 0.69 | 0.014 |
| Rain harvest tank | 1.46 | 0.89 – 2.41 | 0.136 |
| Elevated | 2.54 | 1.36 – 4.76 | 0.004 |
| Flat | 1.92 | 1.45 – 2.54 | <0.001 |
Sobj_cox_fit <- survfit(Sobj_cox)
plot(Sobj_cox_fit)
3.1.2 Cox proporational hazard model, all covariates
Sobj_cox2 <- coxph(Sobj ~ Type +Topo+`_2_2_a_Caretaker` + as.factor(`_3_1_d_primaryuse`) + as.factor(`_2_1_responsible`) +`_4_4_feecollectionsystem`+`_4_17_planOandM` +`_4_18_enoughfundsrepair` +`_4_19_personrepair`+`_4_20_persontrained` + `_5_0_latrineaccess` + `_5_1_a_facilitator`, data=db2)## Warning in fitter(X, Y, strats, offset, init, control, weights = weights, :
## Loglik converged before variable 12 ; beta may be infinite.We should try different types of Cox models later.1
3.1.2.1 Model results, nicer table
stargazer(Sobj_cox2, type = 'html',ci=TRUE, ci.level=0.95, report = 'vc*')| Dependent variable: | |
| Sobj | |
| TypeBorehole mech pump | -0.075 |
| TypeDug well open | -0.848** |
| TypeDug well handpump | -0.136 |
| TypeSpring protected | 0.327 |
| TypeSpring unprotected | -1.444* |
| TypeRain harvest tank | |
| TopoElevated | 0.978** |
| TopoFlat | 0.744*** |
_2_2_a_Caretaker
|
0.409* |
as.factor(_3_1_d_primaryuse)1
|
0.082 |
as.factor(_3_1_d_primaryuse)2
|
0.195 |
as.factor(_3_1_d_primaryuse)3
|
-14.610 |
as.factor(_3_1_d_primaryuse)4
|
0.152 |
as.factor(_2_1_responsible)2
|
0.625 |
as.factor(_2_1_responsible)3
|
-0.314 |
as.factor(_2_1_responsible)4
|
1.420 |
as.factor(_2_1_responsible)5
|
|
as.factor(_2_1_responsible)6
|
0.070 |
as.factor(_2_1_responsible)7
|
-0.352 |
as.factor(_2_1_responsible)8
|
-0.189 |
as.factor(_2_1_responsible)9
|
0.103 |
as.factor(_2_1_responsible)10
|
2.410** |
_4_4_feecollectionsystem
|
-0.365 |
_4_17_planOandM
|
-0.268 |
_4_18_enoughfundsrepair
|
0.131 |
_4_19_personrepair
|
0.466** |
_4_20_persontrained
|
-0.114 |
_5_0_latrineaccess
|
0.011 |
_5_1_a_facilitator
|
-0.162 |
| Observations | 219 |
| R2 | 0.299 |
| Max. Possible R2 | 0.999 |
| Log Likelihood | -737.325 |
| Wald Test | 61.700*** (df = 27) |
| LR Test | 77.735*** (df = 27) |
| Score (Logrank) Test | 73.632*** (df = 27) |
| Note: | p<0.1; p<0.05; p<0.01 |
# stargazer(Sobj_cox2, type = 'text',ci=TRUE, ci.level=0.95, report = 'vc*')3.1.2.2 Same model results, ugly table
summary(Sobj_cox2)## Call:
## coxph(formula = Sobj ~ Type + Topo + `_2_2_a_Caretaker` + as.factor(`_3_1_d_primaryuse`) +
## as.factor(`_2_1_responsible`) + `_4_4_feecollectionsystem` +
## `_4_17_planOandM` + `_4_18_enoughfundsrepair` + `_4_19_personrepair` +
## `_4_20_persontrained` + `_5_0_latrineaccess` + `_5_1_a_facilitator`,
## data = db2)
##
## n= 219, number of events= 159
## (111 observations deleted due to missingness)
##
## coef exp(coef) se(coef) z
## TypeBorehole mech pump -7.491e-02 9.278e-01 3.347e-01 -0.224
## TypeDug well open -8.477e-01 4.284e-01 3.931e-01 -2.157
## TypeDug well handpump -1.363e-01 8.726e-01 2.435e-01 -0.560
## TypeSpring protected 3.268e-01 1.386e+00 3.599e-01 0.908
## TypeSpring unprotected -1.444e+00 2.360e-01 7.784e-01 -1.855
## TypeRain harvest tank NA NA 0.000e+00 NA
## TopoElevated 9.785e-01 2.660e+00 4.176e-01 2.343
## TopoFlat 7.440e-01 2.104e+00 1.910e-01 3.896
## `_2_2_a_Caretaker` 4.086e-01 1.505e+00 2.439e-01 1.675
## as.factor(`_3_1_d_primaryuse`)1 8.242e-02 1.086e+00 1.082e+00 0.076
## as.factor(`_3_1_d_primaryuse`)2 1.951e-01 1.215e+00 2.462e-01 0.793
## as.factor(`_3_1_d_primaryuse`)3 -1.461e+01 4.519e-07 1.827e+03 -0.008
## as.factor(`_3_1_d_primaryuse`)4 1.522e-01 1.164e+00 2.171e-01 0.701
## as.factor(`_2_1_responsible`)2 6.245e-01 1.867e+00 5.149e-01 1.213
## as.factor(`_2_1_responsible`)3 -3.139e-01 7.306e-01 6.797e-01 -0.462
## as.factor(`_2_1_responsible`)4 1.420e+00 4.139e+00 1.145e+00 1.240
## as.factor(`_2_1_responsible`)5 NA NA 0.000e+00 NA
## as.factor(`_2_1_responsible`)6 7.029e-02 1.073e+00 5.021e-01 0.140
## as.factor(`_2_1_responsible`)7 -3.516e-01 7.036e-01 1.169e+00 -0.301
## as.factor(`_2_1_responsible`)8 -1.892e-01 8.276e-01 5.659e-01 -0.334
## as.factor(`_2_1_responsible`)9 1.027e-01 1.108e+00 4.992e-01 0.206
## as.factor(`_2_1_responsible`)10 2.410e+00 1.113e+01 1.161e+00 2.076
## `_4_4_feecollectionsystem` -3.649e-01 6.943e-01 2.417e-01 -1.510
## `_4_17_planOandM` -2.685e-01 7.645e-01 2.101e-01 -1.278
## `_4_18_enoughfundsrepair` 1.306e-01 1.139e+00 2.089e-01 0.625
## `_4_19_personrepair` 4.657e-01 1.593e+00 2.138e-01 2.178
## `_4_20_persontrained` -1.143e-01 8.920e-01 2.877e-01 -0.397
## `_5_0_latrineaccess` 1.115e-02 1.011e+00 7.091e-03 1.572
## `_5_1_a_facilitator` -1.621e-01 8.504e-01 1.975e-01 -0.820
## Pr(>|z|)
## TypeBorehole mech pump 0.8229
## TypeDug well open 0.0310 *
## TypeDug well handpump 0.5756
## TypeSpring protected 0.3639
## TypeSpring unprotected 0.0636 .
## TypeRain harvest tank NA
## TopoElevated 0.0191 *
## TopoFlat 9.79e-05 ***
## `_2_2_a_Caretaker` 0.0939 .
## as.factor(`_3_1_d_primaryuse`)1 0.9393
## as.factor(`_3_1_d_primaryuse`)2 0.4281
## as.factor(`_3_1_d_primaryuse`)3 0.9936
## as.factor(`_3_1_d_primaryuse`)4 0.4832
## as.factor(`_2_1_responsible`)2 0.2252
## as.factor(`_2_1_responsible`)3 0.6442
## as.factor(`_2_1_responsible`)4 0.2149
## as.factor(`_2_1_responsible`)5 NA
## as.factor(`_2_1_responsible`)6 0.8887
## as.factor(`_2_1_responsible`)7 0.7636
## as.factor(`_2_1_responsible`)8 0.7382
## as.factor(`_2_1_responsible`)9 0.8370
## as.factor(`_2_1_responsible`)10 0.0379 *
## `_4_4_feecollectionsystem` 0.1310
## `_4_17_planOandM` 0.2013
## `_4_18_enoughfundsrepair` 0.5319
## `_4_19_personrepair` 0.0294 *
## `_4_20_persontrained` 0.6913
## `_5_0_latrineaccess` 0.1159
## `_5_1_a_facilitator` 0.4120
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## TypeBorehole mech pump 9.278e-01 1.078e+00 0.48149 1.7879
## TypeDug well open 4.284e-01 2.334e+00 0.19826 0.9256
## TypeDug well handpump 8.726e-01 1.146e+00 0.54145 1.4062
## TypeSpring protected 1.386e+00 7.213e-01 0.68479 2.8071
## TypeSpring unprotected 2.360e-01 4.237e+00 0.05134 1.0852
## TypeRain harvest tank NA NA NA NA
## TopoElevated 2.660e+00 3.759e-01 1.17359 6.0307
## TopoFlat 2.104e+00 4.752e-01 1.44729 3.0597
## `_2_2_a_Caretaker` 1.505e+00 6.646e-01 0.93290 2.4268
## as.factor(`_3_1_d_primaryuse`)1 1.086e+00 9.209e-01 0.13022 9.0555
## as.factor(`_3_1_d_primaryuse`)2 1.215e+00 8.227e-01 0.75018 1.9693
## as.factor(`_3_1_d_primaryuse`)3 4.519e-07 2.213e+06 0.00000 Inf
## as.factor(`_3_1_d_primaryuse`)4 1.164e+00 8.588e-01 0.76084 1.7821
## as.factor(`_2_1_responsible`)2 1.867e+00 5.355e-01 0.68070 5.1229
## as.factor(`_2_1_responsible`)3 7.306e-01 1.369e+00 0.19281 2.7682
## as.factor(`_2_1_responsible`)4 4.139e+00 2.416e-01 0.43841 39.0789
## as.factor(`_2_1_responsible`)5 NA NA NA NA
## as.factor(`_2_1_responsible`)6 1.073e+00 9.321e-01 0.40098 2.8703
## as.factor(`_2_1_responsible`)7 7.036e-01 1.421e+00 0.07121 6.9521
## as.factor(`_2_1_responsible`)8 8.276e-01 1.208e+00 0.27298 2.5093
## as.factor(`_2_1_responsible`)9 1.108e+00 9.024e-01 0.41660 2.9478
## as.factor(`_2_1_responsible`)10 1.113e+01 8.986e-02 1.14438 108.2210
## `_4_4_feecollectionsystem` 6.943e-01 1.440e+00 0.43234 1.1149
## `_4_17_planOandM` 7.645e-01 1.308e+00 0.50647 1.1541
## `_4_18_enoughfundsrepair` 1.139e+00 8.776e-01 0.75669 1.7159
## `_4_19_personrepair` 1.593e+00 6.277e-01 1.04778 2.4222
## `_4_20_persontrained` 8.920e-01 1.121e+00 0.50757 1.5677
## `_5_0_latrineaccess` 1.011e+00 9.889e-01 0.99725 1.0254
## `_5_1_a_facilitator` 8.504e-01 1.176e+00 0.57737 1.2525
##
## Concordance= 0.711 (se = 0.022 )
## Rsquare= 0.299 (max possible= 0.999 )
## Likelihood ratio test= 77.73 on 27 df, p=8e-07
## Wald test = 61.7 on 27 df, p=2e-04
## Score (logrank) test = 73.63 on 27 df, p=3e-063.2 Interpretation and commentary
In a typical survival analysis (on, for example, the relationship between amount of medication taken and time until death), hazard ratio (HR) is associated with survival in this way: 2
- HR = 1: No effect
- HR < 1: Reduction in the hazard [meaning longer survival]
- HR > 1: Increase in hazard [meaning shorter survival]
If HR for medicine = 0.5, this means that more medicine means longer survival (compared to less medicine).
In our water quality data, “longer survival” is equivalent to more contamination, so we have to interpret as follows: HR coefficients between 0 and 1 mean positive predicted relationship between the variable and TTC level (more of the variable means more TTC). And HR coefficients above 1 indicate a negative relationship between the variable and TTC level (less of the variable means more TTC).
Interpretation reminders for our data:
- HR < 1: Lower hazard, increased survival, more contamination
- HR > 1: Higher hazard, decreased survival, less contamination
For example, the HR for Type Spring unprotected is 0.1689. That means that the risk of “quicker death” / less contamination is predicted to change by a factor of 0.17 (a reduction3) for unprotected springs, compared to the reference group (boreholes with handpumps). Since the risk of less contamination is predicted to be reduced, this means that *contamination is predicted to be higher in unprotected springs compared to boreholes with handpumps.
This model also predicts that Elevated and Flat topographies are associated with lower levels of contamination than the reference group topography (low lying).
4 Examples of similar applications of survival analysis
- This is more recent and supports the use of Kaplan–Meier method over many others (in some cases), specifically mentioning water quality data: Helsel, Dennis. Much Ado About Next to Nothing: Incorporating Nondetects in Science. The Annals of Occupational Hygiene, Volume 54, Issue 3, April 2010, Pages 257–262, https://doi.org/10.1093/annhyg/mep092. https://academic.oup.com/annweh/article/54/3/257/223531.
- https://onlinelibrary.wiley.com/doi/pdf/10.1111/j.1752-1688.1997.tb03536.x
- https://pubs.er.usgs.gov/publication/70030209
- https://conference.ifas.ufl.edu/GEER2008/Presentation_PDFs/Wednesday/Royal%20Palm%20I-II/1620%20J%20Castro.pdf
- https://link.springer.com/article/10.1007/BF00394476
- Dekker: https://www.sciencedirect.com/science/article/pii/S0085253815534704
- https://www.epa.gov/sites/production/files/2016-05/documents/tech_notes_10_jun2014_r.pdf (suggests MLE over KM)
See: Terry Therneau. Coxph and loglik converged before variable X. https://stat.ethz.ch/pipermail/r-help/2008-September/174201.html.↩
Springate, David. Survival analysis. March 2014. https://rpubs.com/daspringate/survival.↩
Additional guidance on interpretation came from: Cox Proportional-Hazards Model. Statistical tools for high-throughput data analysis. www.sthda.com/english/wiki/cox-proportional-hazards-model.↩