I am using R survival package, survdiff function. I wonder how to get the p value from the return value.
> diff = survdiff(Surv(Time, Censored) ~ Treatment+Gender, data = dat)
> diff
Call:
survdiff(formula = Surv(Time, Censored) ~ Treatment + Gender,
data = dat)
N Observed Expected (O-E)^2/E (O-E)^2/V
Treatment=Control, Gender=M 2 1 1.65 0.255876 0.360905
Treatment=Control, Gender=F 7 3 2.72 0.027970 0.046119
Treatment=IND, Gender=M 5 2 2.03 0.000365 0.000519
Treatment=IND, Gender=F 6 2 1.60 0.100494 0.139041
Chisq= 0.5 on 3 degrees of freedom, p= 0.924
I want to get the p value 0.924 using some function. Thanks.
The code in the function print.survdiff
that displays that values is:
cat("\n Chisq=", format(round(x$chisq, 1)), " on", df,
"degrees of freedom, p=", format(signif(1 - pchisq(x$chisq,
df), digits)), "\n")
The code leading up to it:
if (is.matrix(x$obs)) {
otmp <- apply(x$obs, 1, sum)
etmp <- apply(x$exp, 1, sum)
} else {
otmp <- x$obs
etmp <- x$exp
}
df <- (sum(1 * (etmp > 0))) - 1
And 'digits' is set to 3 in the argument list, so using the example on the surv.diff
help page:
x <- survdiff(Surv(time, status) ~ pat.karno + strata(inst), data=lung)
cat( "p=", format(signif(1 - pchisq(x$chisq,
df), digits)) )
#p= 0.00326
Addressing the comment: In the example the second code block reduces to:
df <- with(x, (sum(1 * (apply(x$exp, 1, sum) > 0))) - 1 )
> df
[1] 7
With the glance()
function from the broom
package, it is very easy to get the p.value.
diff = survdiff(Surv(Time, Censored) ~ Treatment+Gender, data = dat)
broom::glance(diff)$p.value
The p value is not stored in the survdiff class, so it must be calculated on the fly at the time of output. To reproduce the p value one could use the chisq distribution function: "pchisq"
diff = survdiff(Surv(Time, Censored) ~ Treatment+Gender, data = dat)
pchisq(diff$chisq, length(diff$n)-1, lower.tail = FALSE)
If you love us? You can donate to us via Paypal or buy me a coffee so we can maintain and grow! Thank you!
Donate Us With