Calculate hazard from digitized survival curves

user_906592 0
Viewed 65

I have digitized published survival curves and got the following (Survival,time) points

df<-structure(list(Survival = c(1, 1, 1, 1, 1, 1, 0.99, 0.99, 0.99, 
0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 
0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.99, 0.98, 0.98, 0.98, 0.98, 
0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 
0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 0.98, 
0.98, 0.98, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 
0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 
0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 
0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 0.97, 
0.97, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 
0.96, 0.96, 0.96, 0.96, 0.96, 0.96, 0.95, 0.95, 0.95, 0.95, 0.95, 
0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 
0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 
0.95, 0.95, 0.95, 0.95, 0.95, 0.95, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 0.94, 
0.94, 0.94, 0.94, 0.94, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 0.93, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 
0.92, 0.92, 0.92, 0.92, 0.92, 0.92, 0.91, 0.91, 0.91, 0.91, 0.91, 
0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 
0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 
0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 
0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.91, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 0.9, 
0.9, 0.9, 0.9, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 
0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 0.89, 
0.89, 0.89, 0.89, 0.89, 0.89), time = c(0, 0.34, 0.55, 0.82, 
0.87, 0.87, 0.87, 0.87, 0.87, 1.02, 1.02, 1.24, 1.45, 1.61, 1.74, 
1.83, 1.86, 1.87, 1.87, 1.87, 1.88, 1.96, 2.16, 2.37, 2.75, 3.01, 
3.48, 4.27, 4.5, 4.68, 4.89, 5.19, 5.58, 6.43, 6.75, 7.11, 7.55, 
7.99, 8.41, 8.82, 9.26, 9.72, 10.09, 10.41, 10.68, 10.99, 11.29, 
11.58, 11.85, 12.2, 12.47, 12.69, 13.03, 13.31, 13.72, 14.18, 
14.48, 14.79, 15.26, 15.83, 16.08, 16.55, 16.92, 17.74, 18.54, 
19.24, 19.54, 19.92, 21.14, 21.51, 21.89, 22.98, 23.29, 23.55, 
23.78, 24.1, 24.38, 24.79, 25.09, 25.42, 25.78, 27.09, 27.66, 
28.18, 28.57, 28.97, 29.36, 29.76, 30.15, 30.54, 30.94, 31.33, 
31.73, 32.12, 32.64, 33.19, 33.56, 33.8, 34.05, 34.22, 34.5, 
34.79, 35.1, 35.61, 36.05, 36.47, 36.87, 37.26, 37.65, 38.05, 
38.44, 38.84, 39.23, 39.63, 40.02, 40.42, 40.81, 41.21, 41.6, 
41.99, 42.51, 43.1, 43.41, 43.98, 44.91, 45.55, 46.53, 47.44, 
47.83, 48.17, 49.03, 49.65, 50.45, 50.84, 51.22, 51.61, 52.4, 
52.85, 53.24, 53.62, 54.62, 55.11, 55.46, 56, 56.61, 57.39, 57.75, 
58.15, 58.47, 58.77, 59.01, 59.57, 59.94, 60.44, 60.95, 61.35, 
61.74, 62.14, 62.53, 62.93, 63.32, 63.71, 64.11, 64.62, 65.3, 
65.61, 65.96, 66.21, 66.71, 67.08, 67.6, 68.04, 68.45, 68.85, 
69.24, 69.64, 70.03, 70.43, 70.82, 71.22, 71.61, 72.01, 72.4, 
72.79, 73.19, 73.58, 73.98, 74.37, 74.79, 75.23, 75.68, 76.13, 
76.41, 76.65, 76.88, 77.12, 77.52, 77.88, 78.3, 78.8, 79.22, 
79.64, 79.95, 80.94, 81.32, 81.66, 82.5, 82.86, 83.24, 83.58, 
84.01, 84.27, 84.63, 84.94, 85.28, 86.25, 86.59, 86.91, 88.39, 
88.73, 89.3, 89.59, 90.12, 90.39, 90.67, 91.18, 91.7, 92.14, 
92.55, 92.95, 93.34, 93.74, 94.13, 94.53, 94.92, 95.31, 95.71, 
96.1, 96.5, 96.89, 97.29, 97.68, 98.08, 98.47, 98.87, 99.26, 
99.66, 100.05, 100.44, 100.84, 101.23, 101.63, 102.02, 102.42, 
102.81, 103.21, 103.6, 104, 104.39, 104.78, 105.18, 105.57, 105.97, 
106.36, 106.76, 107.15, 107.55, 107.94, 108.34, 108.73, 109.13, 
109.52, 109.91, 110.31, 110.7, 111.12, 111.56, 112.27, 112.64, 
113.1, 113.59, 114.1, 114.46, 114.8, 115.24, 115.55, 116, 116.34, 
117.03, 117.85, 118.24, 118.58, 118.89, 120.23, 120.55, 120.85, 
121.92, 122.4, 122.89, 123.33, 123.75, 124.14, 124.54, 124.93, 
125.33, 125.72, 126.11, 126.51, 126.9, 127.3, 127.69, 128.09, 
128.48, 128.88, 129.27, 129.67, 130.06, 130.45, 130.85, 131.24, 
131.64, 132.03, 132.43, 132.82, 133.22, 133.61, 134.01, 134.4, 
134.8, 135.19, 135.58, 135.98, 136.37, 136.77, 137.16, 137.56, 
137.95, 138.35, 138.74, 139.14, 139.53, 139.92, 140.34, 140.78, 
141.49, 141.86, 142.29, 142.73, 143.3, 143.81, 144.27, 144.73, 
145.23, 145.66, 146.01, 146.33, 146.92, 147.74, 148.26, 148.66, 
149.21, 149.93, 150.33, 150.87, 151.38, 151.78, 152.17, 152.57, 
152.96, 153.36, 153.75, 154.14, 154.54, 154.93, 155.33, 155.72, 
156.12, 156.51, 156.91, 157.3, 157.81, 158.49, 158.82, 159.33, 
159.8, 160.19, 160.74, 161.25, 161.76, 162.44, 162.76, 163.1, 
163.35, 163.98, 164.38, 164.9, 165.25, 165.57, 166.11, 167.09, 
167.48, 167.9, 168.34, 168.75, 169.15, 169.54, 169.94, 170.35, 
170.79, 171.51, 171.88, 172.31, 172.74, 173.32, 173.82, 174.26, 
174.68, 175.07, 175.47, 175.86, 176.26, 176.65, 177.04, 177.44, 
177.83, 178.23, 178.62, 179.02, 179.41, 179.81, 180.2, 180.6, 
180.99, 181.38, 181.78, 182.17, 182.57, 182.96, 183.36, 183.75, 
184.15, 184.54, 184.94, 185.33, 185.73, 186.12, 186.51, 186.91, 
187.3, 187.7, 188.09, 188.49, 188.88, 189.28, 189.67, 190.07, 
190.46, 190.85, 191.25, 191.64, 192.04, 192.43, 192.83, 193.22, 
193.62, 194.01, 194.41, 194.8, 195.2, 195.59, 195.98, 196.38, 
196.77, 197.17, 197.56, 197.96, 198.37, 198.81, 199.62, 199.91, 
200.33, 200.85, 201.34, 201.84, 202.28, 202.7, 203.09, 203.48, 
203.88, 204.27, 204.67, 205.06, 205.46, 205.85, 206.25, 206.64, 
207.04, 207.43, 207.83, 208.22, 208.61, 209.01, 209.4, 209.8, 
210.19, 210.59, 210.98, 211.38, 211.77, 212.17, 212.56, 212.95, 
213.37, 213.81, 214.33, 214.73, 215.63, 216.28, 216.79, 217.3, 
217.69, 218.09, 218.48, 218.88, 219.27, 219.67, 220.06, 220.46, 
220.85, 221.24, 221.64, 222.03, 222.43, 222.82, 223.22, 223.61, 
224.01, 224.4, 224.8, 225.19, 225.59, 225.98, 226.37, 226.77, 
227.16, 227.56, 227.95, 228.35, 228.74, 229.14, 229.53, 229.93, 
230.32, 230.71, 231.11, 231.5, 231.9, 232.29, 232.81, 233.34, 
233.77, 234.34, 235.23, 235.73, 236.24, 236.64, 237.03, 237.43, 
237.82, 238.22, 238.61, 239.4, 239.77, 240.07, 241.47, 242.05, 
242.56, 242.95, 243.35, 243.74, 244.14, 244.53, 244.93, 245.32, 
245.72, 246.11, 246.51, 246.9, 247.29, 247.69, 248.08, 248.48, 
248.87, 249.27, 249.66, 250.06, 250.45, 250.85, 251.24, 251.63, 
252.03, 252.42, 253.4, 253.99, 254.7, 255.15, 255.56, 256.09, 
256.54, 256.87, 257.19, 258.14, 258.48, 258.8, 259.18, 260.18, 
260.56, 261.05, 261.49, 261.9, 262.3, 262.69, 263.09, 263.48, 
264.14, 264.98, 265.49, 265.87, 266.22, 267.66, 268.15, 268.6, 
269.01, 269.41, 269.8, 270.19, 270.59, 270.98, 271.4, 271.84, 
272.36, 272.76, 273.66, 274.31, 274.82, 275.33, 275.72, 276.12, 
276.63, 277.32, 277.67, 277.91, 278.28, 278.54, 278.84, 279.16, 
279.56, 280, 280.73, 281.48, 282.04, 282.61, 283.31, 284.14, 
284.32, 284.44, 284.64, 284.88, 285.15, 285.45, 285.79, 286.17, 
286.72, 287.15, 287.48, 287.82, 288.35, 289.08, 290.04, 290.32, 
290.82, 291.32, 291.7, 292.2, 292.71, 293.42, 293.94, 294.56, 
295.29, 295.61, 296.27, 297.13, 297.45, 297.79, 298.22, 298.54, 
298.88, 299.41, 299.82, 300.22, 300.61, 301.01, 301.4, 301.8, 
302.19, 302.58, 302.98, 303.37, 303.77, 304.16, 304.56, 305.27, 
306.03, 306.39, 306.86, 307.2, 307.62, 307.95, 308.49, 308.9, 
309.3, 309.69, 310.08, 310.48, 310.87, 311.27, 311.66, 312.06, 
312.45, 312.85, 313.24, 313.64, 314.03, 314.43, 314.82, 315.21, 
315.61, 316, 316.4, 316.79, 317.19, 317.58, 317.98, 318.37, 318.77, 
319.16, 319.55, 319.95, 320.34, 320.74, 321.13, 321.53, 321.92, 
322.32, 323.03, 324.13, 324.49, 325.03, 325.77, 326.1, 326.57, 
326.9, 327.43, 327.85, 328.24, 328.64, 329.03, 329.43, 329.82, 
330.22, 330.61, 331, 331.4, 331.79, 332.19, 332.58, 332.98, 333.37, 
333.79, 334.23, 334.7, 335.03, 335.55, 335.89, 336.36, 336.71, 
337.02, 337.47, 338, 338.51, 338.9, 339.3, 339.69, 340.09, 340.48, 
340.88, 341.27, 341.67, 342.06, 342.45, 342.85, 343.24, 343.64, 
344.03, 344.43, 344.82, 345.22, 345.61, 346.01, 346.4, 346.8, 
347.19, 347.58, 347.98, 348.49, 348.88, 349.24, 349.55, 349.89, 
350.32, 350.65, 350.98, 351.51, 352.25, 353.19, 353.49, 353.83, 
354.27, 354.6, 354.93, 355.46, 355.88, 356.27, 356.67, 357.06, 
357.46, 357.85, 358.25, 358.64, 359.3, 360.16, 360.78, 361.24, 
361.79, 362.27, 362.83, 363.19, 363.79)), class = "data.frame", row.names = c(NA, 
-848L))

From which i get

plot(df$time, df$Survival, ylim=c(0,1), type="s", lwd=2, bty="n", xlim=c(0,400))

P

How can I calculate the smoothed instantaneous hazard function from this dataframe? I dont have access to number at risk.

hazard survival-analysis r
1 Answers

The theoretical instantaneous hazard is the time derivative of the cumulative Hazard function which is just 1-Survival. I disagree that the "number at risk" should be necessary in this instance, because you have already created an empirical Survival function. Looking at the Survival curve you would expect a correct estimate to have an early peak and then taper off. So see whether this estimator suffices:

#Pseudo-code
IHF <- { ( S(x+1) - S(x) ) /S(x) }/ # the S(x) in denominator replaces  n-at-risk
            (t(x+1)-t(x)
# R-code
IHF <- function( S, t) { (-diff(S)/head(S, -1))/diff(t) }
plot( head(df$time,-1) , IHF(df$Survival, df$time) )

enter image description here

Well wasn't entirely satisfactory. You only have 11 instances where there is a decrement, so I guess I can see why you want a smoothed estimator. Lets see if loess can suffice. Nope. Couldn't get a sufficient bandwidth. The muhaz package is designed for this. I'm not sure if it's the best for it, but it's the one I'm most familiar with. It designed for censored data but the default for its censoring status is to assume no cases are censored, so yo can just provide times and it will take the survival into account:

fit2 <- muhaz(head(df$time,-1)[ as.logical(diff(df$Survival)) ], maxtime=360)
plot(fit2)

enter image description here

I wasn't entirely happy with that result, being concerned about the late peak. So I tried anotehr package:

install.packages("bshazard"
library(bshazard
times=head(df$time,-1)[ as.logical(diff(df$Survival)) ]
 times
#[1]   0.87   3.48  13.72  33.80  75.68  90.39 145.66 200.85 272.76 295.29 352.25
  fit<-bshazard( Surv(times,rep(1, length(times))) ~ 1, data=df)

enter image description here

I line that one a bit better. (Which is a bit ironic since I'm the maintainer of pkg:muhaz, but I'm just a caretaker.)

user_18556770
Related