Here is the exact code I have written which does the standard vs nt1 and
standard vs nt2 and also gives me the hazard ratio for nt1 vs nt2.

with <- read.table("allwiths.txt",header=TRUE)
fix(arm)
    function (data)
    {
    dummy <- rep(0,2437)
    for(i in 1:2437){
    if(data$Arm[i]=="CBZ")
    dummy[i] <- i
    }
    return(data[-dummy,])
    }
fix(x1)
    function (data)
    {
    x1 <- rep(0,716)
    for(i in 1:716){
    if(data$Treat[i]=="TPM")
    x1[i] <- 0
    if(data$Treat[i]=="VPS")
    x1[i] <- 0
    if(data$Treat[i]=="LTG")
    x1[i] <- 1
    }
    return(x1)
    }
fix(x2)
    function (data)
    {
    x2 <- rep(0,716)
    for(i in 1:716){
    if(data$Treat[i]=="TPM")
    x2[i] <- 1
    if(data$Treat[i]=="VPS")
    x2[i] <- 0
    if(data$Treat[i]=="LTG")
    x2[i] <- 0
    }
    return(x2)
    }
fit1 = coxph(Surv(Withtime,Wcens)~x1(armb),data=armb,method="breslow")
exp(fit1$coefficients)
exp(confint(fit1))
fit2 = coxph(Surv(Withtime,Wcens)~x2(armb),data=armb,method="breslow")
exp(fit2$coefficients)
exp(confint(fit2))
exp(fit2$coefficients)/exp(fit1$coefficients)

From that, how do I get the necessary variance-covaraince matrix.

Sorry if I appear dense.  It really isn't my intention.

Laura


On Wed, Aug 27, 2008 at 10:36 AM, Peter Dalgaard
<P.Dalgaard@biostat.ku.dk>wrote:

> Laura Bonnett wrote:
> > Hi all,
> >
> > Sorry to ask again but I'm still not sure how to get the full
> > variance-covariance matrix.  Peter suggested a three-level treatment
> > factor.  However, I thought that the censoring variable could only take
> > values 0 or 1 so how do you programme such a factor.
> >
> Well, not to put it too diplomatically: If you can confuse the treatment
> factor and the censoring indicator, there might be more wrong with your
> analysis than we originally assumed, so how about showing us a bit more
> of what you actually did?
>
> > Alternatively, is there another way to produce the required covariance?
> >
> > Thank you,
> >
> > Laura
> >
> > On Tue, Aug 26, 2008 at 11:37 AM, Laura Bonnett
> > <l.j.bonnett@googlemail.com>wrote:
> >
> >
> >> The standard treatment is the same in both comparison.
> >>
> >> How do you do a three-level treatment factor?
> >> I thought you had to have a censoring indicator which took values 0 or 1
> >> not 1, 2 or 3?
> >>
> >> Thanks,
> >>
> >> Laura
> >>
> >> On Tue, Aug 26, 2008 at 11:05 AM, Peter Dalgaard <
> P.Dalgaard@biostat.ku.dk
> >>
> >>> wrote:
> >>>
> >>> Laura Bonnett wrote:
> >>>
> >>>> Dear R help forum,
> >>>>
> >>>> I am using the function 'coxph' to obtain hazard ratios for the
> >>>>
> >>> comparison
> >>>
> >>>> of a standard treatment to new treatments.  This is easily obtained by
> >>>> fitting the relevant model and then calling exp(coef(fit1)) say.
> >>>>
> >>>> I now want to obtain the hazard ratio for the comparison of two
> >>>>
> >>> non-standard
> >>>
> >>>> treatments.
> >>>> >From a statistical point of view, this can be achieved by dividing
> the
> >>>> exponentiated coefficients of 2 comparisions. E.g. to compared new
> >>>>
> >>> treatment
> >>>
> >>>> 1 (nt1) to new treatment 2 (nt2) we can fit 2 models:
> >>>> fit1 = standard treatment vs nt1
> >>>> fit2 = standard treatment vs nt2.
> >>>> The required hazard ratio is therefore exp(coef(fit1))/exp(coef(fit2))
> >>>>
> >>>> In order to obtain an associated confidence interval for this I
> require
> >>>>
> >>> the
> >>>
> >>>> covariance of this comparison.  I know that R gives the
> >>>>
> >>> variance-covariance
> >>>
> >>>> matrix by the command 'fit$var'.  However, this only gives the
> >>>>
> >>> covariance
> >>>
> >>>> matrix for non standard drugs and not the full covariance matrix.
> >>>>
> >>>> Can anyone tell me how to obtain the full covariance matrix?
> >>>>
> >>>>
> >>>>
> >>> What kind of data do you have? Is the "standard treatment group" the
> >>> same in both comparisons?  If so, why not just have a three-level
> >>> treatment factor and compare nt1 to nt2 directly. If the control groups
> >>> are completely separate, then the covariance between fits made on
> >>> independent data is of course zero.
> >>>
> >>>
> >>>> Thank you,
> >>>>
> >>>> Laura
> >>>>
> >>>>       [[alternative HTML version deleted]]
> >>>>
> >>>> ______________________________________________
> >>>> R-help@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.
> >>>>
> >>>>
> >>> --
> >>>   O__  ---- Peter Dalgaard             �ster Farimagsgade 5, Entr.B
> >>>  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
> >>>  (*) \(*) -- University of Copenhagen   Denmark      Ph:  (+45)
> 35327918
> >>> ~~~~~~~~~~ - (p.dalgaard@biostat.ku.dk)              FAX: (+45)
> 35327907
> >>>
> >>>
> >>>
> >>>
> >
> >       [[alternative HTML version deleted]]
> >
> >
> > ------------------------------------------------------------------------
> >
> > ______________________________________________
> > R-help@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.
> >
>
>
> --
>   O__  ---- Peter Dalgaard             Øster Farimagsgade 5, Entr.B
>  c/ /'_ --- Dept. of Biostatistics     PO Box 2099, 1014 Cph. K
>  (*) \(*) -- University of Copenhagen   Denmark      Ph:  (+45) 35327918
> ~~~~~~~~~~ - (p.dalgaard@biostat.ku.dk)              FAX: (+45) 35327907
>
>
>

	[[alternative HTML version deleted]]

