# Analysis of killdeer data using the proportional hazerds model
# 2019-05-01 CJS Initial code
# This is the killdeer data that ships with RMark
library(AICcmodavg)
library(ggplot2)
## Registered S3 methods overwritten by 'ggplot2':
## method from
## [.quosures rlang
## c.quosures rlang
## print.quosures rlang
library(plyr)
library(readxl)
library(survival)
library(survminer)
## Loading required package: ggpubr
## Loading required package: magrittr
##
## Attaching package: 'ggpubr'
## The following object is masked from 'package:plyr':
##
## mutate
source(file.path("..","..","logistic-exposure-model.R")) # need to expand the nest data
# The dataframe must contain the following fields with the following names
#
# FirstFound: day the nest was first found
# LastPresent: last day that a chick was present in the nest
# LastChecked: last day the nest was checked
# Fate: fate of the nest; 0=hatch an
# Freq: number of nests with this data
# AgeDay1 - age at day 1 so that age of best at each day can be imputed
#
# In this example, the multiple visits to a nest have been collapsed
# to a single record for each nest.
# In more complex examples, you may have multple records per nest
# as shown in the mallard example.
#
killdata <- readxl::read_excel(file.path("..","Killdeer.xlsx"),
sheet="killdeer-age")
head(killdata)
## # A tibble: 6 x 7
## id FirstFound LastPresent LastChecked Fate Freq AgeDay1
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 /*A*/ 1 9 9 0 1 0
## 2 /*B*/ 5 5 9 1 1 -2
## 3 /*C*/ 5 40 40 0 1 -3
## 4 /*D*/ 9 32 32 0 1 -4
## 5 /*E*/ 7 8 8 0 1 -4
## 6 /*F*/ 3 15 15 0 1 1
# There are no covariates except nest age
killdata2 <- expand.nest.data.ph(killdata)
head(killdata2)
## id FirstFound LastPresent LastChecked Fate Freq AgeDay1 Start End
## 1 /*A*/ 1 9 9 0 1 0 1 2
## 2 /*A*/ 1 9 9 0 1 0 2 3
## 3 /*A*/ 1 9 9 0 1 0 3 4
## 4 /*A*/ 1 9 9 0 1 0 4 5
## 5 /*A*/ 1 9 9 0 1 0 5 6
## 6 /*A*/ 1 9 9 0 1 0 6 7
## Fail NestAge Surv
## 1 0 0 (1,2+]
## 2 0 1 (2,3+]
## 3 0 2 (3,4+]
## 4 0 3 (4,5+]
## 5 0 4 (5,6+]
## 6 0 5 (6,7+]
killdata2
## id FirstFound LastPresent LastChecked Fate Freq AgeDay1 Start End
## 1 /*A*/ 1 9 9 0 1 0 1 2
## 2 /*A*/ 1 9 9 0 1 0 2 3
## 3 /*A*/ 1 9 9 0 1 0 3 4
## 4 /*A*/ 1 9 9 0 1 0 4 5
## 5 /*A*/ 1 9 9 0 1 0 5 6
## 6 /*A*/ 1 9 9 0 1 0 6 7
## 7 /*A*/ 1 9 9 0 1 0 7 8
## 8 /*A*/ 1 9 9 0 1 0 8 9
## 9 /*B*/ 5 5 9 1 1 -2 5 9
## 10 /*C*/ 5 40 40 0 1 -3 5 6
## 11 /*C*/ 5 40 40 0 1 -3 6 7
## 12 /*C*/ 5 40 40 0 1 -3 7 8
## 13 /*C*/ 5 40 40 0 1 -3 8 9
## 14 /*C*/ 5 40 40 0 1 -3 9 10
## 15 /*C*/ 5 40 40 0 1 -3 10 11
## 16 /*C*/ 5 40 40 0 1 -3 11 12
## 17 /*C*/ 5 40 40 0 1 -3 12 13
## 18 /*C*/ 5 40 40 0 1 -3 13 14
## 19 /*C*/ 5 40 40 0 1 -3 14 15
## 20 /*C*/ 5 40 40 0 1 -3 15 16
## 21 /*C*/ 5 40 40 0 1 -3 16 17
## 22 /*C*/ 5 40 40 0 1 -3 17 18
## 23 /*C*/ 5 40 40 0 1 -3 18 19
## 24 /*C*/ 5 40 40 0 1 -3 19 20
## 25 /*C*/ 5 40 40 0 1 -3 20 21
## 26 /*C*/ 5 40 40 0 1 -3 21 22
## 27 /*C*/ 5 40 40 0 1 -3 22 23
## 28 /*C*/ 5 40 40 0 1 -3 23 24
## 29 /*C*/ 5 40 40 0 1 -3 24 25
## 30 /*C*/ 5 40 40 0 1 -3 25 26
## 31 /*C*/ 5 40 40 0 1 -3 26 27
## 32 /*C*/ 5 40 40 0 1 -3 27 28
## 33 /*C*/ 5 40 40 0 1 -3 28 29
## 34 /*C*/ 5 40 40 0 1 -3 29 30
## 35 /*C*/ 5 40 40 0 1 -3 30 31
## 36 /*C*/ 5 40 40 0 1 -3 31 32
## 37 /*C*/ 5 40 40 0 1 -3 32 33
## 38 /*C*/ 5 40 40 0 1 -3 33 34
## 39 /*C*/ 5 40 40 0 1 -3 34 35
## 40 /*C*/ 5 40 40 0 1 -3 35 36
## 41 /*C*/ 5 40 40 0 1 -3 36 37
## 42 /*C*/ 5 40 40 0 1 -3 37 38
## 43 /*C*/ 5 40 40 0 1 -3 38 39
## 44 /*C*/ 5 40 40 0 1 -3 39 40
## 45 /*D*/ 9 32 32 0 1 -4 9 10
## 46 /*D*/ 9 32 32 0 1 -4 10 11
## 47 /*D*/ 9 32 32 0 1 -4 11 12
## 48 /*D*/ 9 32 32 0 1 -4 12 13
## 49 /*D*/ 9 32 32 0 1 -4 13 14
## 50 /*D*/ 9 32 32 0 1 -4 14 15
## 51 /*D*/ 9 32 32 0 1 -4 15 16
## 52 /*D*/ 9 32 32 0 1 -4 16 17
## 53 /*D*/ 9 32 32 0 1 -4 17 18
## 54 /*D*/ 9 32 32 0 1 -4 18 19
## 55 /*D*/ 9 32 32 0 1 -4 19 20
## 56 /*D*/ 9 32 32 0 1 -4 20 21
## 57 /*D*/ 9 32 32 0 1 -4 21 22
## 58 /*D*/ 9 32 32 0 1 -4 22 23
## 59 /*D*/ 9 32 32 0 1 -4 23 24
## 60 /*D*/ 9 32 32 0 1 -4 24 25
## 61 /*D*/ 9 32 32 0 1 -4 25 26
## 62 /*D*/ 9 32 32 0 1 -4 26 27
## 63 /*D*/ 9 32 32 0 1 -4 27 28
## 64 /*D*/ 9 32 32 0 1 -4 28 29
## 65 /*D*/ 9 32 32 0 1 -4 29 30
## 66 /*D*/ 9 32 32 0 1 -4 30 31
## 67 /*D*/ 9 32 32 0 1 -4 31 32
## 68 /*E*/ 7 8 8 0 1 -4 7 8
## 69 /*F*/ 3 15 15 0 1 1 3 4
## 70 /*F*/ 3 15 15 0 1 1 4 5
## 71 /*F*/ 3 15 15 0 1 1 5 6
## 72 /*F*/ 3 15 15 0 1 1 6 7
## 73 /*F*/ 3 15 15 0 1 1 7 8
## 74 /*F*/ 3 15 15 0 1 1 8 9
## 75 /*F*/ 3 15 15 0 1 1 9 10
## 76 /*F*/ 3 15 15 0 1 1 10 11
## 77 /*F*/ 3 15 15 0 1 1 11 12
## 78 /*F*/ 3 15 15 0 1 1 12 13
## 79 /*F*/ 3 15 15 0 1 1 13 14
## 80 /*F*/ 3 15 15 0 1 1 14 15
## 81 /*G*/ 8 32 32 0 1 -7 8 9
## 82 /*G*/ 8 32 32 0 1 -7 9 10
## 83 /*G*/ 8 32 32 0 1 -7 10 11
## 84 /*G*/ 8 32 32 0 1 -7 11 12
## 85 /*G*/ 8 32 32 0 1 -7 12 13
## 86 /*G*/ 8 32 32 0 1 -7 13 14
## 87 /*G*/ 8 32 32 0 1 -7 14 15
## 88 /*G*/ 8 32 32 0 1 -7 15 16
## 89 /*G*/ 8 32 32 0 1 -7 16 17
## 90 /*G*/ 8 32 32 0 1 -7 17 18
## 91 /*G*/ 8 32 32 0 1 -7 18 19
## 92 /*G*/ 8 32 32 0 1 -7 19 20
## 93 /*G*/ 8 32 32 0 1 -7 20 21
## 94 /*G*/ 8 32 32 0 1 -7 21 22
## 95 /*G*/ 8 32 32 0 1 -7 22 23
## 96 /*G*/ 8 32 32 0 1 -7 23 24
## 97 /*G*/ 8 32 32 0 1 -7 24 25
## 98 /*G*/ 8 32 32 0 1 -7 25 26
## 99 /*G*/ 8 32 32 0 1 -7 26 27
## 100 /*G*/ 8 32 32 0 1 -7 27 28
## 101 /*G*/ 8 32 32 0 1 -7 28 29
## 102 /*G*/ 8 32 32 0 1 -7 29 30
## 103 /*G*/ 8 32 32 0 1 -7 30 31
## 104 /*G*/ 8 32 32 0 1 -7 31 32
## 105 /*H*/ 14 14 16 1 1 -10 14 16
## 106 /*I*/ 8 14 14 0 1 -7 8 9
## 107 /*I*/ 8 14 14 0 1 -7 9 10
## 108 /*I*/ 8 14 14 0 1 -7 10 11
## 109 /*I*/ 8 14 14 0 1 -7 11 12
## 110 /*I*/ 8 14 14 0 1 -7 12 13
## 111 /*I*/ 8 14 14 0 1 -7 13 14
## 112 /*J*/ 13 14 14 0 1 -12 13 14
## 113 /*K*/ 14 33 33 0 1 -13 14 15
## 114 /*K*/ 14 33 33 0 1 -13 15 16
## 115 /*K*/ 14 33 33 0 1 -13 16 17
## 116 /*K*/ 14 33 33 0 1 -13 17 18
## 117 /*K*/ 14 33 33 0 1 -13 18 19
## 118 /*K*/ 14 33 33 0 1 -13 19 20
## 119 /*K*/ 14 33 33 0 1 -13 20 21
## 120 /*K*/ 14 33 33 0 1 -13 21 22
## 121 /*K*/ 14 33 33 0 1 -13 22 23
## 122 /*K*/ 14 33 33 0 1 -13 23 24
## 123 /*K*/ 14 33 33 0 1 -13 24 25
## 124 /*K*/ 14 33 33 0 1 -13 25 26
## 125 /*K*/ 14 33 33 0 1 -13 26 27
## 126 /*K*/ 14 33 33 0 1 -13 27 28
## 127 /*K*/ 14 33 33 0 1 -13 28 29
## 128 /*K*/ 14 33 33 0 1 -13 29 30
## 129 /*K*/ 14 33 33 0 1 -13 30 31
## 130 /*K*/ 14 33 33 0 1 -13 31 32
## 131 /*K*/ 14 33 33 0 1 -13 32 33
## 132 /*L*/ 15 37 37 0 1 -10 15 16
## 133 /*L*/ 15 37 37 0 1 -10 16 17
## 134 /*L*/ 15 37 37 0 1 -10 17 18
## 135 /*L*/ 15 37 37 0 1 -10 18 19
## 136 /*L*/ 15 37 37 0 1 -10 19 20
## 137 /*L*/ 15 37 37 0 1 -10 20 21
## 138 /*L*/ 15 37 37 0 1 -10 21 22
## 139 /*L*/ 15 37 37 0 1 -10 22 23
## 140 /*L*/ 15 37 37 0 1 -10 23 24
## 141 /*L*/ 15 37 37 0 1 -10 24 25
## 142 /*L*/ 15 37 37 0 1 -10 25 26
## 143 /*L*/ 15 37 37 0 1 -10 26 27
## 144 /*L*/ 15 37 37 0 1 -10 27 28
## 145 /*L*/ 15 37 37 0 1 -10 28 29
## 146 /*L*/ 15 37 37 0 1 -10 29 30
## 147 /*L*/ 15 37 37 0 1 -10 30 31
## 148 /*L*/ 15 37 37 0 1 -10 31 32
## 149 /*L*/ 15 37 37 0 1 -10 32 33
## 150 /*L*/ 15 37 37 0 1 -10 33 34
## 151 /*L*/ 15 37 37 0 1 -10 34 35
## 152 /*L*/ 15 37 37 0 1 -10 35 36
## 153 /*L*/ 15 37 37 0 1 -10 36 37
## 154 /*M*/ 16 37 40 1 1 -11 16 17
## 155 /*M*/ 16 37 40 1 1 -11 17 18
## 156 /*M*/ 16 37 40 1 1 -11 18 19
## 157 /*M*/ 16 37 40 1 1 -11 19 20
## 158 /*M*/ 16 37 40 1 1 -11 20 21
## 159 /*M*/ 16 37 40 1 1 -11 21 22
## 160 /*M*/ 16 37 40 1 1 -11 22 23
## 161 /*M*/ 16 37 40 1 1 -11 23 24
## 162 /*M*/ 16 37 40 1 1 -11 24 25
## 163 /*M*/ 16 37 40 1 1 -11 25 26
## 164 /*M*/ 16 37 40 1 1 -11 26 27
## 165 /*M*/ 16 37 40 1 1 -11 27 28
## 166 /*M*/ 16 37 40 1 1 -11 28 29
## 167 /*M*/ 16 37 40 1 1 -11 29 30
## 168 /*M*/ 16 37 40 1 1 -11 30 31
## 169 /*M*/ 16 37 40 1 1 -11 31 32
## 170 /*M*/ 16 37 40 1 1 -11 32 33
## 171 /*M*/ 16 37 40 1 1 -11 33 34
## 172 /*M*/ 16 37 40 1 1 -11 34 35
## 173 /*M*/ 16 37 40 1 1 -11 35 36
## 174 /*M*/ 16 37 40 1 1 -11 36 37
## 175 /*M*/ 16 37 40 1 1 -11 37 40
## 176 /*N*/ 16 28 32 1 1 -11 16 17
## 177 /*N*/ 16 28 32 1 1 -11 17 18
## 178 /*N*/ 16 28 32 1 1 -11 18 19
## 179 /*N*/ 16 28 32 1 1 -11 19 20
## 180 /*N*/ 16 28 32 1 1 -11 20 21
## 181 /*N*/ 16 28 32 1 1 -11 21 22
## 182 /*N*/ 16 28 32 1 1 -11 22 23
## 183 /*N*/ 16 28 32 1 1 -11 23 24
## 184 /*N*/ 16 28 32 1 1 -11 24 25
## 185 /*N*/ 16 28 32 1 1 -11 25 26
## 186 /*N*/ 16 28 32 1 1 -11 26 27
## 187 /*N*/ 16 28 32 1 1 -11 27 28
## 188 /*N*/ 16 28 32 1 1 -11 28 32
## 189 /*O*/ 16 17 17 0 1 -12 16 17
## 190 /*P*/ 21 28 33 1 1 -15 21 22
## 191 /*P*/ 21 28 33 1 1 -15 22 23
## 192 /*P*/ 21 28 33 1 1 -15 23 24
## 193 /*P*/ 21 28 33 1 1 -15 24 25
## 194 /*P*/ 21 28 33 1 1 -15 25 26
## 195 /*P*/ 21 28 33 1 1 -15 26 27
## 196 /*P*/ 21 28 33 1 1 -15 27 28
## 197 /*P*/ 21 28 33 1 1 -15 28 33
## 198 /*Q*/ 23 33 34 1 1 -17 23 24
## 199 /*Q*/ 23 33 34 1 1 -17 24 25
## 200 /*Q*/ 23 33 34 1 1 -17 25 26
## 201 /*Q*/ 23 33 34 1 1 -17 26 27
## 202 /*Q*/ 23 33 34 1 1 -17 27 28
## 203 /*Q*/ 23 33 34 1 1 -17 28 29
## 204 /*Q*/ 23 33 34 1 1 -17 29 30
## 205 /*Q*/ 23 33 34 1 1 -17 30 31
## 206 /*Q*/ 23 33 34 1 1 -17 31 32
## 207 /*Q*/ 23 33 34 1 1 -17 32 33
## 208 /*Q*/ 23 33 34 1 1 -17 33 34
## 209 /*R*/ 27 29 29 0 1 -23 27 28
## 210 /*R*/ 27 29 29 0 1 -23 28 29
## Fail NestAge Surv
## 1 0 0 ( 1, 2+]
## 2 0 1 ( 2, 3+]
## 3 0 2 ( 3, 4+]
## 4 0 3 ( 4, 5+]
## 5 0 4 ( 5, 6+]
## 6 0 5 ( 6, 7+]
## 7 0 6 ( 7, 8+]
## 8 0 7 ( 8, 9+]
## 9 1 2 ( 5, 9]
## 10 0 1 ( 5, 6+]
## 11 0 2 ( 6, 7+]
## 12 0 3 ( 7, 8+]
## 13 0 4 ( 8, 9+]
## 14 0 5 ( 9,10+]
## 15 0 6 (10,11+]
## 16 0 7 (11,12+]
## 17 0 8 (12,13+]
## 18 0 9 (13,14+]
## 19 0 10 (14,15+]
## 20 0 11 (15,16+]
## 21 0 12 (16,17+]
## 22 0 13 (17,18+]
## 23 0 14 (18,19+]
## 24 0 15 (19,20+]
## 25 0 16 (20,21+]
## 26 0 17 (21,22+]
## 27 0 18 (22,23+]
## 28 0 19 (23,24+]
## 29 0 20 (24,25+]
## 30 0 21 (25,26+]
## 31 0 22 (26,27+]
## 32 0 23 (27,28+]
## 33 0 24 (28,29+]
## 34 0 25 (29,30+]
## 35 0 26 (30,31+]
## 36 0 27 (31,32+]
## 37 0 28 (32,33+]
## 38 0 29 (33,34+]
## 39 0 30 (34,35+]
## 40 0 31 (35,36+]
## 41 0 32 (36,37+]
## 42 0 33 (37,38+]
## 43 0 34 (38,39+]
## 44 0 35 (39,40+]
## 45 0 4 ( 9,10+]
## 46 0 5 (10,11+]
## 47 0 6 (11,12+]
## 48 0 7 (12,13+]
## 49 0 8 (13,14+]
## 50 0 9 (14,15+]
## 51 0 10 (15,16+]
## 52 0 11 (16,17+]
## 53 0 12 (17,18+]
## 54 0 13 (18,19+]
## 55 0 14 (19,20+]
## 56 0 15 (20,21+]
## 57 0 16 (21,22+]
## 58 0 17 (22,23+]
## 59 0 18 (23,24+]
## 60 0 19 (24,25+]
## 61 0 20 (25,26+]
## 62 0 21 (26,27+]
## 63 0 22 (27,28+]
## 64 0 23 (28,29+]
## 65 0 24 (29,30+]
## 66 0 25 (30,31+]
## 67 0 26 (31,32+]
## 68 0 2 ( 7, 8+]
## 69 0 3 ( 3, 4+]
## 70 0 4 ( 4, 5+]
## 71 0 5 ( 5, 6+]
## 72 0 6 ( 6, 7+]
## 73 0 7 ( 7, 8+]
## 74 0 8 ( 8, 9+]
## 75 0 9 ( 9,10+]
## 76 0 10 (10,11+]
## 77 0 11 (11,12+]
## 78 0 12 (12,13+]
## 79 0 13 (13,14+]
## 80 0 14 (14,15+]
## 81 0 0 ( 8, 9+]
## 82 0 1 ( 9,10+]
## 83 0 2 (10,11+]
## 84 0 3 (11,12+]
## 85 0 4 (12,13+]
## 86 0 5 (13,14+]
## 87 0 6 (14,15+]
## 88 0 7 (15,16+]
## 89 0 8 (16,17+]
## 90 0 9 (17,18+]
## 91 0 10 (18,19+]
## 92 0 11 (19,20+]
## 93 0 12 (20,21+]
## 94 0 13 (21,22+]
## 95 0 14 (22,23+]
## 96 0 15 (23,24+]
## 97 0 16 (24,25+]
## 98 0 17 (25,26+]
## 99 0 18 (26,27+]
## 100 0 19 (27,28+]
## 101 0 20 (28,29+]
## 102 0 21 (29,30+]
## 103 0 22 (30,31+]
## 104 0 23 (31,32+]
## 105 1 3 (14,16]
## 106 0 0 ( 8, 9+]
## 107 0 1 ( 9,10+]
## 108 0 2 (10,11+]
## 109 0 3 (11,12+]
## 110 0 4 (12,13+]
## 111 0 5 (13,14+]
## 112 0 0 (13,14+]
## 113 0 0 (14,15+]
## 114 0 1 (15,16+]
## 115 0 2 (16,17+]
## 116 0 3 (17,18+]
## 117 0 4 (18,19+]
## 118 0 5 (19,20+]
## 119 0 6 (20,21+]
## 120 0 7 (21,22+]
## 121 0 8 (22,23+]
## 122 0 9 (23,24+]
## 123 0 10 (24,25+]
## 124 0 11 (25,26+]
## 125 0 12 (26,27+]
## 126 0 13 (27,28+]
## 127 0 14 (28,29+]
## 128 0 15 (29,30+]
## 129 0 16 (30,31+]
## 130 0 17 (31,32+]
## 131 0 18 (32,33+]
## 132 0 4 (15,16+]
## 133 0 5 (16,17+]
## 134 0 6 (17,18+]
## 135 0 7 (18,19+]
## 136 0 8 (19,20+]
## 137 0 9 (20,21+]
## 138 0 10 (21,22+]
## 139 0 11 (22,23+]
## 140 0 12 (23,24+]
## 141 0 13 (24,25+]
## 142 0 14 (25,26+]
## 143 0 15 (26,27+]
## 144 0 16 (27,28+]
## 145 0 17 (28,29+]
## 146 0 18 (29,30+]
## 147 0 19 (30,31+]
## 148 0 20 (31,32+]
## 149 0 21 (32,33+]
## 150 0 22 (33,34+]
## 151 0 23 (34,35+]
## 152 0 24 (35,36+]
## 153 0 25 (36,37+]
## 154 0 4 (16,17+]
## 155 0 5 (17,18+]
## 156 0 6 (18,19+]
## 157 0 7 (19,20+]
## 158 0 8 (20,21+]
## 159 0 9 (21,22+]
## 160 0 10 (22,23+]
## 161 0 11 (23,24+]
## 162 0 12 (24,25+]
## 163 0 13 (25,26+]
## 164 0 14 (26,27+]
## 165 0 15 (27,28+]
## 166 0 16 (28,29+]
## 167 0 17 (29,30+]
## 168 0 18 (30,31+]
## 169 0 19 (31,32+]
## 170 0 20 (32,33+]
## 171 0 21 (33,34+]
## 172 0 22 (34,35+]
## 173 0 23 (35,36+]
## 174 0 24 (36,37+]
## 175 1 25 (37,40]
## 176 0 4 (16,17+]
## 177 0 5 (17,18+]
## 178 0 6 (18,19+]
## 179 0 7 (19,20+]
## 180 0 8 (20,21+]
## 181 0 9 (21,22+]
## 182 0 10 (22,23+]
## 183 0 11 (23,24+]
## 184 0 12 (24,25+]
## 185 0 13 (25,26+]
## 186 0 14 (26,27+]
## 187 0 15 (27,28+]
## 188 1 16 (28,32]
## 189 0 3 (16,17+]
## 190 0 5 (21,22+]
## 191 0 6 (22,23+]
## 192 0 7 (23,24+]
## 193 0 8 (24,25+]
## 194 0 9 (25,26+]
## 195 0 10 (26,27+]
## 196 0 11 (27,28+]
## 197 1 12 (28,33]
## 198 0 5 (23,24+]
## 199 0 6 (24,25+]
## 200 0 7 (25,26+]
## 201 0 8 (26,27+]
## 202 0 9 (27,28+]
## 203 0 10 (28,29+]
## 204 0 11 (29,30+]
## 205 0 12 (30,31+]
## 206 0 13 (31,32+]
## 207 0 14 (32,33+]
## 208 1 15 (33,34]
## 209 0 3 (27,28+]
## 210 0 4 (28,29+]
# Caution -- Age and Time are very highly correlated!
ggplot(data=killdata2, aes(x=Start, y=NestAge))+
geom_point()
fit.Sdot.ph <- coxph(Surv ~1, data=killdata2)
summary(fit.Sdot.ph)
## Call: coxph(formula = Surv ~ 1, data = killdata2)
##
## Null model
## log likelihood= -9.651945
## n= 210
fit.SNestAge.ph <- coxph(Surv~NestAge, data=killdata2)
summary(fit.SNestAge.ph)
## Call:
## coxph(formula = Surv ~ NestAge, data = killdata2)
##
## n= 210, number of events= 6
##
## coef exp(coef) se(coef) z Pr(>|z|)
## NestAge -0.3359 0.7147 0.1728 -1.944 0.0519 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## exp(coef) exp(-coef) lower .95 upper .95
## NestAge 0.7147 1.399 0.5094 1.003
##
## Concordance= 0.815 (se = 0.084 )
## Likelihood ratio test= 6.95 on 1 df, p=0.008
## Wald test = 3.78 on 1 df, p=0.05
## Score (logrank) test = 5.53 on 1 df, p=0.02
cox.zph(fit.SNestAge.ph)
## rho chisq p
## NestAge -0.493 0.428 0.513
ggcoxzph(cox.zph(fit.SNestAge.ph))
png(file=file.path("..","..","..","..","MyStuff","Images","killdeer-ph-age-testprop.png"), h=4, w=6, units="in", res=300)
ggcoxzph(cox.zph(fit.SNestAge.ph))
dev.off()
## quartz_off_screen
## 2
# baseline cumulative hazard
cumhaz <- basehaz(fit.SNestAge.ph)
# estimate change in cumulative hazard and plot
cumhaz$deltaHaz <- c(NA,diff(cumhaz$hazard))
ggplot(data=cumhaz, aes(x=time, y=deltaHaz))+
ggtitle("Estimated baseline hazard function over time")+
geom_point()+
geom_smooth(se=FALSE)
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## Warning: Removed 1 rows containing non-finite values (stat_smooth).
## Warning: Removed 1 rows containing missing values (geom_point).
# estimated hazard curves at three ages
pred.data <- expand.grid( NestAge=seq(min(killdata2$NestAge),max(killdata2$NestAge), length.out=3))
pred.survival <- survfit(fit.SNestAge.ph, newdata=pred.data, se.fit=TRUE)
plot(pred.survival)
AICcmodavg::aictab(list(fit.Sdot.ph, fit.SNestAge.ph))
## Warning in aictab.AICcoxph(list(fit.Sdot.ph, fit.SNestAge.ph)):
## Model names have been supplied automatically in the table
##
## Model selection based on AICc:
##
## K AICc Delta_AICc AICcWt Cum.Wt LL
## Mod2 1 14.38 0.00 0.92 0.92 -6.18
## Mod1 0 19.30 4.93 0.08 1.00 -9.65
# Now to fit the logistic exposure model
killdata3 <- expand.nest.data(killdata)
fit.Sdot <- glm(Survive~1,
family=binomial(link=logexp(killdata3$Exposure)),
data=killdata3)
fit.SNestAge <- glm(Survive~NestAge,
family=binomial(link=logexp(killdata3$Exposure)),
data=killdata3)
AICcmodavg::aictab(list(fit.Sdot, fit.SNestAge))
## Warning in aictab.AICglm.lm(list(fit.Sdot, fit.SNestAge)):
## Model names have been supplied automatically in the table
##
## Model selection based on AICc:
##
## K AICc Delta_AICc AICcWt Cum.Wt LL
## Mod1 1 44.53 0.00 0.71 0.71 -21.26
## Mod2 2 46.31 1.79 0.29 1.00 -21.13