[Rd] Possible bug in summary.survfit - 'scale' argument ignored?
Marc Schwartz
marc_schwartz at me.com
Fri Apr 3 18:12:55 CEST 2009
On Apr 1, 2009, at 7:08 AM, Thomas Lumley wrote:
>
> I've sent a fixed version 2.35-4 to CRAN. It turned out to be a
> fairly simple change.
>
> -thomas
Thomas,
Apologies for the delay in my reply. I am on vacation and will only
have sporadic e-mail access thru mid-next week.
I noted the fix that you made and it is indeed pretty simple. I do
have one comment/question however.
It would seem that the nature of the interaction between the 'scale'
and 'times' arguments has now changed with this update. Presumably as
a result of the transform taking place at the beginning of the old
version as opposed to being at the end of the new function.
For example, using the prior version of survival included with R 2.8.1
(2.34-1):
# Get the summary output in months (eg. days / 30.44)
# Note that 'times' is in the 'scale' transformed based time interval
(0:21 months)
> summary( survfit( Surv(futime, fustat)~1, data=ovarian), scale =
30.44, times = 0:21)
Call: survfit(formula = Surv(futime, fustat) ~ 1, data = ovarian)
time n.risk n.event survival std.err lower 95% CI upper 95% CI
0 26 0 1.000 0.0000 1.000 1.000
1 26 0 1.000 0.0000 1.000 1.000
2 25 1 0.962 0.0377 0.890 1.000
3 25 0 0.962 0.0377 0.890 1.000
4 24 1 0.923 0.0523 0.826 1.000
5 24 0 0.923 0.0523 0.826 1.000
6 23 1 0.885 0.0627 0.770 1.000
7 23 0 0.885 0.0627 0.770 1.000
8 23 0 0.885 0.0627 0.770 1.000
9 22 1 0.846 0.0708 0.718 0.997
10 22 0 0.846 0.0708 0.718 0.997
11 21 1 0.808 0.0773 0.670 0.974
12 19 2 0.731 0.0870 0.579 0.923
13 18 0 0.731 0.0870 0.579 0.923
14 17 0 0.731 0.0870 0.579 0.923
15 15 1 0.688 0.0919 0.529 0.894
16 12 2 0.596 0.0999 0.429 0.828
17 12 0 0.596 0.0999 0.429 0.828
18 12 0 0.596 0.0999 0.429 0.828
19 11 1 0.546 0.1032 0.377 0.791
20 11 0 0.546 0.1032 0.377 0.791
21 10 1 0.497 0.1051 0.328 0.752
However, in the new version:
> summary( survfit( Surv(futime, fustat)~1, data=ovarian), scale =
30.44, times = 0:21)
Call: survfit(formula = Surv(futime, fustat) ~ 1, data = ovarian)
time n.risk n.event survival std.err lower 95% CI upper 95% CI
0.0000 26 0 1 0 1 1
0.0329 26 0 1 0 1 1
0.0657 26 0 1 0 1 1
0.0986 26 0 1 0 1 1
0.1314 26 0 1 0 1 1
0.1643 26 0 1 0 1 1
0.1971 26 0 1 0 1 1
0.2300 26 0 1 0 1 1
0.2628 26 0 1 0 1 1
0.2957 26 0 1 0 1 1
0.3285 26 0 1 0 1 1
0.3614 26 0 1 0 1 1
0.3942 26 0 1 0 1 1
0.4271 26 0 1 0 1 1
0.4599 26 0 1 0 1 1
0.4928 26 0 1 0 1 1
0.5256 26 0 1 0 1 1
0.5585 26 0 1 0 1 1
0.5913 26 0 1 0 1 1
0.6242 26 0 1 0 1 1
0.6570 26 0 1 0 1 1
0.6899 26 0 1 0 1 1
It would appear that 'times' are 0:21 days, but as fractions of a month:
> (0:21) / 30.44
[1] 0.00000000 0.03285151 0.06570302 0.09855453 0.13140604 0.16425756
[7] 0.19710907 0.22996058 0.26281209 0.29566360 0.32851511 0.36136662
[13] 0.39421813 0.42706965 0.45992116 0.49277267 0.52562418 0.55847569
[19] 0.59132720 0.62417871 0.65703022 0.68988173
To get the same output as in the prior version:
> summary( survfit( Surv(futime, fustat)~1, data=ovarian), scale =
30.44, times = (0:21) * 30.44)
Call: survfit(formula = Surv(futime, fustat) ~ 1, data = ovarian)
time n.risk n.event survival std.err lower 95% CI upper 95% CI
0 26 0 1.000 0.0000 1.000 1.000
1 26 0 1.000 0.0000 1.000 1.000
2 25 1 0.962 0.0377 0.890 1.000
3 25 0 0.962 0.0377 0.890 1.000
4 24 1 0.923 0.0523 0.826 1.000
5 24 0 0.923 0.0523 0.826 1.000
6 23 1 0.885 0.0627 0.770 1.000
7 23 0 0.885 0.0627 0.770 1.000
8 23 0 0.885 0.0627 0.770 1.000
9 22 1 0.846 0.0708 0.718 0.997
10 22 0 0.846 0.0708 0.718 0.997
11 21 1 0.808 0.0773 0.670 0.974
12 19 2 0.731 0.0870 0.579 0.923
13 18 0 0.731 0.0870 0.579 0.923
14 17 0 0.731 0.0870 0.579 0.923
15 15 1 0.688 0.0919 0.529 0.894
16 12 2 0.596 0.0999 0.429 0.828
17 12 0 0.596 0.0999 0.429 0.828
18 12 0 0.596 0.0999 0.429 0.828
19 11 1 0.546 0.1032 0.377 0.791
20 11 0 0.546 0.1032 0.377 0.791
21 10 1 0.497 0.1051 0.328 0.752
Truth be told, the help for summary.survfit() is not explicit on this
point, but presumably, there is much code out in the wild based upon
the prior behavior.
I hope this is helpful.
Regards,
Marc
More information about the R-devel
mailing list