[R] print and coef Methods for survreg Differ
biii m@iii@g oii de@@ey@ws
biii m@iii@g oii de@@ey@ws
Tue Feb 23 15:45:53 CET 2021
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]]
More information about the R-help
mailing list