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

Steve Lianoglou mailinglist.honeypot at gmail.com
Wed Jul 6 00:22:38 CEST 2011


Quick note:

On Tue, Jul 5, 2011 at 6:16 PM, Bert Gunter <gunter.berton at gene.com> wrote:
> Yes, it can be done using basic plot commands.
>
> But if you really want to get fancy and plot "grouped" graphs, I
> strongly recommend you look into R's packages -- ggplot or trellis.

Attempting to clear out any confusion before it sets in: I'm pretty
sure Bert meant "lattice" instead of "trellis".

-steve

> Both have excellent documentation and companion books and  were built
> for this sort of thing. The (considerable) learning curve will be
> worth the effort.
>
> Cheers,
> Bert
>
> On Tue, Jul 5, 2011 at 3:08 PM, Trey Batey <ekt.batey at gmail.com> 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.
>>
>> --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.
>>
>
>
>
> --
> "Men by nature long to get on to the ultimate truths, and will often
> be impatient with elementary studies or fight shy of them. If it were
> possible to reach the ultimate truths without the elementary studies
> usually prefixed to them, these would not be preparatory studies but
> superfluous diversions."
>
> -- Maimonides (1135-1204)
>
> Bert Gunter
> Genentech Nonclinical Biostatistics
> 467-7374
> http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm
>
> ______________________________________________
> 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.
>



-- 
Steve Lianoglou
Graduate Student: Computational Systems Biology
 | Memorial Sloan-Kettering Cancer Center
 | Weill Medical College of Cornell University
Contact Info: http://cbio.mskcc.org/~lianos/contact



More information about the R-help mailing list