# This file reproduces Figure S2 rm(list=ls(all=TRUE)) load("/data/seer/pops.RData") # this loads in pops pops1=subset(pops,popyear<1985) pops2=subset(pops,(popyear>=1985)&(popyear<1997)) pops3=subset(pops,popyear>=1997) getSexPY<-function(pops) { attach(pops) pym=NULL;pyf=NULL for (i in 0:18) { pym[i+1]=sum(population[(popsex==1)&(popage==i)]) pyf[i+1]=sum(population[(popsex==2)&(popage==i)]) } detach(pops) list(m=pym,f=pyf) } (py1=getSexPY(pops1)) (py2=getSexPY(pops2)) (py3=getSexPY(pops3)) load("/data/SEER/lymyleuk.RData") # this loads in DF #if (k==8) { canc=c("NHL", "Hodgkins", "MML" ,"ALL"," CLL", "AML", "CML") # code=list(c(2000,2008),c(2015,2016),c(2030,2030),c(2040,2040),c(2041,2041),c(2050,2050),c(2051,2051))} #Indx=(DF$ICD9==2041)&(DF$numPrim==1);tit="CLL" #Indx=(DF$ICD9==2051)&(DF$numPrim==1);tit="CML" #including CMML, not as clear #Indx=(DF$morph==9866)&(DF$numPrim==1);tit="APL" Indx=(DF$morph==9827)&(DF$numPrim==1);tit="ATL" #Indx=(DF$ICD9==2050)&(DF$numPrim==1);tit="AML" #Indx=(DF$ICD9<2009)&(DF$ICD9>=2000)&(DF$numPrim==1);tit="NHLL" #Indx=(DF$ICD9==2040)&(DF$numPrim==1); tit="ALL" #perhaps males are consistent with CML, but noisy #Indx=(DF$ICD9==2030)&(DF$numPrim==1); tit="MML" ndat=DF[Indx,] ndat1=subset(ndat,yDiag<1985) ndat2=subset(ndat,(yDiag>=1985)&(yDiag<1997)) ndat3=subset(ndat,yDiag>=1997) age=c(0.5,3,seq(7.5,87.5,5)) m1=with(ndat1,hist(ageRC[sex==1],breaks=c(seq(-.5,17.5,1),100),plot=FALSE)$counts) f1=with(ndat1,hist(ageRC[sex==2],breaks=c(seq(-.5,17.5,1),100),plot=FALSE)$counts) m2=with(ndat2,hist(ageRC[sex==1],breaks=c(seq(-.5,17.5,1),100),plot=FALSE)$counts) f2=with(ndat2,hist(ageRC[sex==2],breaks=c(seq(-.5,17.5,1),100),plot=FALSE)$counts) m3=with(ndat3,hist(ageRC[sex==1],breaks=c(seq(-.5,17.5,1),100),plot=FALSE)$counts) f3=with(ndat3,hist(ageRC[sex==2],breaks=c(seq(-.5,17.5,1),100),plot=FALSE)$counts) graphics.off() windows(height=5,width=12) par(mfrow=c(1,2),mar=c(4.7,6.2,0.8,1),lwd=3,cex.lab=2,cex.axis=2,cex.main=1.3) #par(mfcol=c(1,2),mar=c(4.5,4.1,.2,.2),cex.lab=1.4) plot(age,m1/py1$m,log="y",type='b',xlab="Age",pch=1,col='blue',ylab="",yaxt="n",cex=1.3) lines(age,m2/py2$m,type='b',col='red',pch=2,cex=1.3) lines(age,m3/py3$m,type='b',col='black',pch=3,cex=1.3) axis(side=2,las=1, at=c(1e-7,1e-6,1e-5,1e-4),labels=expression(.1,1,10,10^2)) #mtext(expression("Cases per 10^6 Person-Years"),side=2,line=3.5,cex=2) mtext(expression(paste("Cases per ",10^6," PY")),side=2,line=3.7,cex=2) mtext(paste(tit,"Males"),side=3,line=-1.8,cex=2,adj=0.05) legend("bottomright",c("1973-1984","1985-1996","1997-2008"), col=c("blue","red","black"),pch=1:3,bty="n",cex=1.5) plot(age,f1/py1$f,log="y",type='b',xlab="Age",pch=1,cex=1.3,col='blue',ylab="",yaxt="n" ) lines(age,f2/py2$f,type='b',col='red',pch=2,cex=1.3) lines(age,f3/py3$f,type='b',col='black',pch=3,cex=1.3) mtext(paste(tit,"Females"),side=3,line=-1.8,cex=2,adj=0.05) axis(side=2,las=1, at=c(1e-7,1e-6,1e-5,1e-4),labels=expression(.1,1,10,10^2)) mtext(expression(paste("Cases per ",10^6," PY")),side=2,line=3.5,cex=2)