[R] print and coef Methods for survreg Differ
Jeff Newmiller
jdnewm|| @end|ng |rom dcn@d@v|@@c@@u@
Tue Feb 23 18:24:27 CET 2021
You can't. But in order to be able to predict for states that _were_ in the training data, the coefficient cannot be NA.
On February 23, 2021 8:40:43 AM PST, bill using denney.ws wrote:
>How should you be able to make a prediction (using this type of model)
>from a state where there is no data such as treatment="C" in my
>example?
>
>-----Original Message-----
>From: Jeff Newmiller <jdnewmil using dcn.davis.ca.us>
>Sent: Tuesday, February 23, 2021 11:10 AM
>To: r-help using r-project.org; bill using denney.ws
>Subject: Re: [R] print and coef Methods for survreg Differ
>
>Model equations do not normally have conditional forms dependent on
>whether specific coefficients are NA or not. If you assign NA to a
>coefficient then you will not be able to predict outputs for input
>cases that you should be able to. Zero allows these expected cases to
>work... NA would prevent any useful prediction output.
>
>On February 23, 2021 6:45:53 AM PST, bill using denney.ws wrote:
>>Hello,
>>
>>
>>
>>I'm working on a survreg model where the full data are subset for
>>modeling individual parts of the data separately. When subsetting,
>the
>>fit variable ("treatment" in the example below) has levels that are
>not
>>in the data.
>> A
>>work-around for this is to drop the levels, but it seems inaccurate to
>
>>have the `coef()` method provide zero as the coefficient for the level
>
>>without data.
>>
>>
>>
>>Why does coef(model) provide zero as the coefficient for treatment
>>instead of NA? Is this a bug?
>>
>>
>>
>>Thanks,
>>
>>
>>
>>Bill
>>
>>
>>
>>``` r
>>
>>library(survival)
>>
>>library(emmeans)
>>
>>
>>
>>my_data <-
>>
>> data.frame(
>>
>> value=c(rep(1, 5), 6:10),
>>
>> treatment=factor(rep(c("A", "B"), each=5), levels=c("A", "B",
>"C"))
>>
>> )
>>
>>my_data$cens <- c(0, 1)[(my_data$value == 1) + 1]
>>
>>
>>
>>model <- survreg(Surv(time=value, event=cens)~treatment, data=my_data)
>>
>>#> Warning in survreg.fit(X, Y, weights, offset, init = init,
>>controlvals =
>>
>>#> control, : Ran out of iterations and did not converge
>>
>>coef(model)
>>
>>#> (Intercept) treatmentB treatmentC
>>
>>#> 0.08588218 2.40341893 0.00000000
>>
>>model$coef
>>
>>#> (Intercept) treatmentB treatmentC
>>
>>#> 0.08588218 2.40341893 NA
>>
>>model$coefficients
>>
>>#> (Intercept) treatmentB treatmentC
>>
>>#> 0.08588218 2.40341893 0.00000000
>>
>>print(model)
>>
>>#> Call:
>>
>>#> survreg(formula = Surv(time = value, event = cens) ~ treatment,
>>
>>#> data = my_data)
>>
>>#>
>>
>>#> Coefficients: (1 not defined because of singularities)
>>
>>#> (Intercept) treatmentB treatmentC
>>
>>#> 0.08588218 2.40341893 NA
>>
>>#>
>>
>>#> Scale= 0.09832254
>>
>>#>
>>
>>#> Loglik(model)= 4.9 Loglik(intercept only)= -15
>>
>>#> Chisq= 39.92 on 2 degrees of freedom, p= 2.15e-09
>>
>>#> n= 10
>>
>>summary(model)
>>
>>#>
>>
>>#> Call:
>>
>>#> survreg(formula = Surv(time = value, event = cens) ~ treatment,
>>
>>#> data = my_data)
>>
>>#> Value Std. Error z p
>>
>>#> (Intercept) 0.0859 0.0681 1.26 0.21
>>
>>#> treatmentB 2.4034 0.2198 10.93 <2e-16
>>
>>#> treatmentC 0.0000 0.0000 NA NA
>>
>>#> Log(scale) -2.3195 0.0000 -Inf <2e-16
>>
>>#>
>>
>>#> Scale= 0.0983
>>
>>#>
>>
>>#> Weibull distribution
>>
>>#> Loglik(model)= 4.9 Loglik(intercept only)= -15
>>
>>#> Chisq= 39.92 on 2 degrees of freedom, p= 2.1e-09
>>
>>#> Number of Newton-Raphson Iterations: 30
>>
>>#> n= 10
>>
>>ref_grid(model)
>>
>>#> Error in ref_grid(model): Something went wrong:
>>
>>#> Non-conformable elements in reference grid.
>>
>>
>>
>>my_data_correct_levels <- my_data
>>
>>my_data_correct_levels$treatment <-
>>droplevels(my_data_correct_levels$treatment)
>>
>>
>>
>>model_correct <- survreg(Surv(time=value, event=cens)~treatment,
>>data=my_data_correct_levels)
>>
>>#> Warning in survreg.fit(X, Y, weights, offset, init = init,
>>controlvals =
>>
>>#> control, : Ran out of iterations and did not converge
>>
>>coef(model_correct)
>>
>>#> (Intercept) treatmentB
>>
>>#> 0.08588218 2.40341893
>>
>>print(model_correct)
>>
>>#> Call:
>>
>>#> survreg(formula = Surv(time = value, event = cens) ~ treatment,
>>
>>#> data = my_data_correct_levels)
>>
>>#>
>>
>>#> Coefficients:
>>
>>#> (Intercept) treatmentB
>>
>>#> 0.08588218 2.40341893
>>
>>#>
>>
>>#> Scale= 0.09832254
>>
>>#>
>>
>>#> Loglik(model)= 4.9 Loglik(intercept only)= -15
>>
>>#> Chisq= 39.92 on 1 degrees of freedom, p= 2.65e-10
>>
>>#> n= 10
>>
>>summary(model_correct)
>>
>>#>
>>
>>#> Call:
>>
>>#> survreg(formula = Surv(time = value, event = cens) ~ treatment,
>>
>>#> data = my_data_correct_levels)
>>
>>#> Value Std. Error z p
>>
>>#> (Intercept) 0.0859 0.0681 1.26 0.21
>>
>>#> treatmentB 2.4034 0.2198 10.93 <2e-16
>>
>>#> Log(scale) -2.3195 0.0000 -Inf <2e-16
>>
>>#>
>>
>>#> Scale= 0.0983
>>
>>#>
>>
>>#> Weibull distribution
>>
>>#> Loglik(model)= 4.9 Loglik(intercept only)= -15
>>
>>#> Chisq= 39.92 on 1 degrees of freedom, p= 2.6e-10
>>
>>#> Number of Newton-Raphson Iterations: 30
>>
>>#> n= 10
>>
>>ref_grid(model_correct)
>>
>>#> 'emmGrid' object with variables:
>>
>>#> treatment = A, B
>>
>>#> Transformation: "log"
>>
>>```
>>
>>
>>
>><sup>Created on 2021-02-23 by the [reprex
>>package](https://reprex.tidyverse.org) (v1.0.0)</sup>
>>
>>
>> [[alternative HTML version deleted]]
>>
>>______________________________________________
>>R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
>>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.
>
>--
>Sent from my phone. Please excuse my brevity.
--
Sent from my phone. Please excuse my brevity.
More information about the R-help
mailing list