[R] plotting survival curves (multiple curves on single graph)

David Winsemius dwinsemius at comcast.net
Wed Jul 6 00:26:48 CEST 2011


On Jul 5, 2011, at 6:24 PM, David Winsemius wrote:

>
> On Jul 5, 2011, at 6:08 PM, Trey Batey wrote:
>
>> Hello.
>>
>> This is a follow-up to a question I posted last week.  With some
>> previous suggestions from the R-help community, I have been able to
>> plot survival (, hazard, and density) curves using published data for
>> Siler hazard parameters from a number of ethnographic populations.
>> Can the function below be modified, perhaps with a "for" statement,  
>> so
>> that multiple curves (different line types---one for each population)
>> are plotted on a single graph for comparison?  Thanks so much.
>>
>
> There are (at least) three methods to plot multiple curves in base  
> plotting:
> -- plot() then lines()
>
> ?lines
>
> --plot(); par(add=TRUE); plot()
>
Er, ... make that par(new=TRUE)

> ?par
>
> # There is also matplot()
> ?matplot
>
> After extracting the sil function to exist on its own, you could try:
>
> matplot(x=t, y=apply(hazanth[ ,3:7], 1, sil)
>
> My first choice would be to make a modified version of your silsurv  
> that uses the lines function rather than plot and then you can just  
> use the lines of code you already have.
>
> -- 
> David
>> --Trey
>>
>> The function and calls below use the data in this Excel file (feel
>> free to access):
>> https://docs.google.com/leaf?id=0B5zZGW2utJN0ZDk1NjA0ZjUtMWU0ZS00ZGQ3LWIxZTUtOWE0NGVmYWMxODJl&hl=en_US
>>
>> ## - plot Siler survival curve
>> ##############################
>> silsurv<-function(a1,b1,a2,a3,b3)
>> {
>>   sil=function(t)
>>     {
>>       h.t<-a1*exp(-b1*t)+a2+a3*exp(b3*t)
>>       S.t<-exp(-a1/b1*(1-exp(-b1*t))-a2*t+a3/b3*(1-exp(b3*t)))
>>       d.t<-S.t*h.t
>>
>>       #return(d.t)
>>       return(S.t)
>>       #return(h.t)
>>     }
>>   t<-seq(0,90,1)
>> plot 
>> (t 
>> ,sil 
>> (t 
>> ),ylim 
>> =c(0,1),type='l',cex.lab=0.8,cex.axis=0.75,ylab='S(t)',xlab='Age
>> (years)')
>> }
>>
>> with 
>> (hazanth 
>> [1,3 
>> : 
>> 7 
>> ],silsurv 
>> (a1 
>> =a1,b1=b1,a2=a2,a3=a3,b3=b3));title(main=hazanth[1,1],cex.main=0.9)
>> # plot for Hadza
>> with 
>> (hazanth 
>> [2,3 
>> : 
>> 7 
>> ],silsurv 
>> (a1 
>> =a1,b1=b1,a2=a2,a3=a3,b3=b3));title(main=hazanth[2,1],cex.main=0.9)
>> # plot for Ache
>> with 
>> (hazanth 
>> [3,3 
>> : 
>> 7 
>> ],silsurv 
>> (a1 
>> =a1,b1=b1,a2=a2,a3=a3,b3=b3));title(main=hazanth[3,1],cex.main=0.9)
>> # plot for Hiwi
>> with 
>> (hazanth 
>> [4,3 
>> : 
>> 7 
>> ],silsurv 
>> (a1 
>> =a1,b1=b1,a2=a2,a3=a3,b3=b3));title(main=hazanth[4,1],cex.main=0.9)
>> # plot for !Kung
>> with 
>> (hazanth 
>> [5,3 
>> : 
>> 7 
>> ],silsurv 
>> (a1 
>> =a1,b1=b1,a2=a2,a3=a3,b3=b3));title(main=hazanth[5,1],cex.main=0.9)
>> # plot for Yanomamo
>> with 
>> (hazanth 
>> [6,3 
>> : 
>> 7 
>> ],silsurv 
>> (a1 
>> =a1,b1=b1,a2=a2,a3=a3,b3=b3));title(main=hazanth[6,1],cex.main=0.9)
>> # plot for Tsimane
>>
>> ###############################
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>
> David Winsemius, MD
> West Hartford, CT
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.

David Winsemius, MD
West Hartford, CT



More information about the R-help mailing list