# for later on... # T4_31146035.txt downloaded from # http://www-dep.iarc.fr/WHOdb/WHOdb.htm # read tab-delimited file deaths=read.table("T4_31146035.txt", skip=2,sep="\t", as.is=rep(TRUE,21)) # avoid levels # look at first 8 lines deaths[1:8,] n.deaths.breast.ca = as.numeric(deaths[6,3:20]) str(n.deaths.breast.ca) mid.age = c(seq(2.5,82.5,5), 90) mean.age.at.br.ca.death =sum(mid.age*n.deaths.breast.ca)/ sum(n.deaths.breast.ca) round(mean.age.at.br.ca.death,1) # look up life expectancy at this age... # or (more refined) get a vector of the 18 life # expectancies closest to the 18 mid-ages # and take a weighted average of them, with weights # given by the n.deaths.breast.ca vector # ................................... # Fig 13.1 # http://www.cancer.ca/ # canada-wide/ # about%20cancer/cancer%20statistics/ x=scan() 19.365 0.978 22.966 2.934 25.736 5.868 27.953 12.714 31.000 21.516 35.986 51.834 37.649 65.526 42.912 120.293 47.345 174.083 52.608 235.697 57.041 277.751 59.811 302.201 61.196 312.958 61.750 318.826 66.736 341.320 72.000 349.144 75.324 358.924 77.264 363.814 81.696 358.924 85.574 354.034 87.236 351.100 89.730 340.342 92.500 326.650 # br.ca.incidence = matrix(x,ncol=2,byrow=TRUE) ; br.ca.incidence x=scan() 21.858 0.978 25.182 0.978 28.784 1.956 31.277 2.934 33.770 5.868 37.095 9.780 41.250 15.648 44.020 20.538 47.068 26.406 49.284 34.230 52.054 44.010 56.209 50.856 59.534 61.614 61.750 70.416 64.243 76.284 66.736 81.174 68.953 89.976 72.000 101.711 73.939 110.513 75.601 119.315 77.541 130.073 78.926 138.875 81.142 156.479 83.358 177.995 84.466 191.687 86.128 208.313 88.345 239.609 89.176 254.279 90.284 276.773 91.392 296.333 92.500 314.914 # br.ca.mortality = matrix(x,ncol=2,byrow=TRUE) ; br.ca.mortality # all cause death rates # http://www.bdlc.umontreal.ca/CHMD/prov/can/can.htm # Canada, Death rates (period 1x1) Last modified: 06-Jul-2010, MPv4 (Feb05) # Year Age Female Male Total x=scan() 2007 0 0.004779 0.005637 0.005219 2007 1 0.000290 0.000329 0.000310 2007 2 0.000188 0.000238 0.000214 2007 3 0.000100 0.000178 0.000140 2007 4 0.000100 0.000111 0.000106 2007 5 0.000142 0.000100 0.000120 2007 6 0.000134 0.000121 0.000127 2007 7 0.000143 0.000129 0.000136 2007 8 0.000095 0.000096 0.000096 2007 9 0.000104 0.000104 0.000104 2007 10 0.000069 0.000091 0.000080 2007 11 0.000081 0.000107 0.000094 2007 12 0.000114 0.000146 0.000131 2007 13 0.000107 0.000172 0.000140 2007 14 0.000181 0.000241 0.000212 2007 15 0.000214 0.000293 0.000254 2007 16 0.000196 0.000381 0.000291 2007 17 0.000258 0.000613 0.000440 2007 18 0.000336 0.000733 0.000540 2007 19 0.000311 0.000894 0.000611 2007 20 0.000295 0.000974 0.000644 2007 21 0.000281 0.000787 0.000541 2007 22 0.000314 0.000791 0.000559 2007 23 0.000291 0.000768 0.000535 2007 24 0.000305 0.000817 0.000566 2007 25 0.000353 0.000776 0.000567 2007 26 0.000307 0.000825 0.000569 2007 27 0.000301 0.000811 0.000558 2007 28 0.000332 0.000728 0.000531 2007 29 0.000330 0.000883 0.000607 2007 30 0.000362 0.000804 0.000583 2007 31 0.000347 0.000817 0.000583 2007 32 0.000512 0.000754 0.000633 2007 33 0.000410 0.000866 0.000639 2007 34 0.000463 0.000972 0.000720 2007 35 0.000472 0.001076 0.000777 2007 36 0.000490 0.001132 0.000814 2007 37 0.000714 0.001100 0.000909 2007 38 0.000580 0.001218 0.000902 2007 39 0.000744 0.001371 0.001060 2007 40 0.000785 0.001353 0.001072 2007 41 0.000936 0.001577 0.001260 2007 42 0.001074 0.001578 0.001328 2007 43 0.001049 0.001794 0.001425 2007 44 0.001153 0.001889 0.001524 2007 45 0.001520 0.002121 0.001822 2007 46 0.001589 0.002291 0.001941 2007 47 0.001792 0.002677 0.002236 2007 48 0.001854 0.002722 0.002288 2007 49 0.002070 0.003024 0.002547 2007 50 0.002238 0.003504 0.002870 2007 51 0.002331 0.003710 0.003017 2007 52 0.002750 0.004358 0.003548 2007 53 0.002882 0.004575 0.003720 2007 54 0.003051 0.005001 0.004015 2007 55 0.003301 0.005488 0.004383 2007 56 0.003748 0.005989 0.004857 2007 57 0.004111 0.006863 0.005471 2007 58 0.004515 0.006962 0.005723 2007 59 0.004691 0.007587 0.006119 2007 60 0.005491 0.009257 0.007346 2007 61 0.005447 0.008904 0.007148 2007 62 0.006001 0.009876 0.007906 2007 63 0.006949 0.011298 0.009085 2007 64 0.007954 0.012343 0.010103 2007 65 0.008344 0.013591 0.010902 2007 66 0.009115 0.014404 0.011680 2007 67 0.010354 0.016330 0.013237 2007 68 0.011641 0.018737 0.015048 2007 69 0.012071 0.020525 0.016106 2007 70 0.013876 0.021915 0.017686 2007 71 0.015365 0.023785 0.019326 2007 72 0.016635 0.026203 0.021113 2007 73 0.017794 0.029660 0.023323 2007 74 0.020296 0.032816 0.026096 2007 75 0.021909 0.036187 0.028458 2007 76 0.024400 0.038978 0.030992 2007 77 0.027681 0.044966 0.035371 2007 78 0.030024 0.048631 0.038161 2007 79 0.034677 0.055322 0.043530 2007 80 0.038098 0.059793 0.047165 2007 81 0.042324 0.064134 0.051186 2007 82 0.049939 0.074766 0.059738 2007 83 0.053699 0.083155 0.064973 2007 84 0.063637 0.091710 0.074076 2007 85 0.072948 0.102935 0.083789 2007 86 0.080629 0.109213 0.090624 2007 87 0.094759 0.130488 0.106741 2007 88 0.098653 0.140145 0.111993 2007 89 0.117613 0.158768 0.130415 2007 90 0.127075 0.171828 0.140460 2007 91 0.150923 0.193049 0.163003 2007 92 0.168475 0.212616 0.180594 2007 93 0.184488 0.231694 0.196765 2007 94 0.209158 0.250106 0.219275 2007 95 0.228212 0.274841 0.239144 2007 96 0.242532 0.289219 0.252815 2007 97 0.270056 0.349879 0.286560 2007 98 0.307209 0.370569 0.319284 2007 99 0.331630 0.415572 0.346670 2007 100 0.359387 0.466732 0.377421 2007 101 0.396426 0.463078 0.406423 2007 102 0.395641 0.468620 0.405474 2007 103 0.483298 0.363482 0.467648 2007 104 0.426042 0.609789 0.449661 2007 105 0.600411 0.541462 0.593322 2007 106 0.581653 0.834299 0.608625 2007 107 0.769546 0.740436 0.766318 # all.cause.mortality = matrix(x,ncol=5,byrow=TRUE)[,3] ; all.cause.mortality ######## # smooth functions, from 25 to 92 incidence = function(age) approx(br.ca.incidence[,1], br.ca.incidence[,2]/100000, age)$y mortality = function(age) approx(br.ca.mortality[,1], br.ca.mortality[,2]/100000, age)$y all.mortality = function(age) approx(0:107,all.cause.mortality, age)$y other.mortality = function(age) all.mortality(age) - mortality(age) plot(25:92, other.mortality(25:92), type="l", log="y" ) points(25:92, incidence(25:92), type="l" , log="y" ) # lifetime risk, with lifetime taken as from 22 to 92 at.risk.this.interval = c() br.ca.cases.this.interval = c() deaths.this.interval = c() candidates = 100 delta.age=1 for (age in seq(22,92,delta.age) ) { at.risk.this.interval = c(at.risk.this.interval, candidates) new.cases = candidates * (1 - exp(-incidence(age+delta.age/2) * delta.age) ) new.deaths = candidates * (1 - exp(-other.mortality(age+delta.age/2) * delta.age) ) candidates = candidates - (new.cases + new.deaths ) br.ca.cases.this.interval = c(br.ca.cases.this.interval, new.cases ) deaths.this.interval = c(deaths.this.interval, new.deaths ) } round(cbind(at.risk.this.interval,br.ca.cases.this.interval, deaths.this.interval ),1) round(cbind(at.risk.this.interval, cumsum(br.ca.cases.this.interval), cumsum(deaths.this.interval) ),3) # other parts to question # Say (adult) life starts at 25 # average longevity from 25 onwards # if mortality rate from breast cancer set to zero, # i.e., if only the other.mortality function were operating.. # get S curves, and E[longevity from 25 to 92] S.without = c(1,exp( - cumsum(other.mortality(25:92)*1) )) plot(24:92,S.without,type="l", col="green", xlim=c(0,95),ylim=c(0,1)) abline(v=24) text(50, .65, toString( round(sum(S.without*1),1) ),col="green",adj=c(1,0)) S.with = c(1,exp( - cumsum(all.mortality(25:92)*1) )) lines(24:92,S.with,col="red") text(50, .50, toString( round(sum(S.with*1),1) ),col="red",adj=c(1,0))