    * Rep_Survival 11 years test of standard errors.sps.
    * Produce file of survival times for teeth filled in period.
    * Adjusted for time dependency.
    * Test standard errors by using random censoring as alternative to weighted adjustment.
    * Written by PSKL on 24/09/02.

    * Work from file of observed restoration lives.

    get file='D:\Longitudinal Data\Rep_Observed Restoration Lives 1991 to 2001.sav'.

    * Restrict to restorations placed before the end of the observation period 31/12/2001.

    select if (doplacn<yrmoda(2001,12,31)-yrmoda(1899,12,31)).

    * Check whether time is in acceptable range, and set reint flag.

    compute reint=0.
    if not (missing(time)) reint=1.
    if (time le 0) reint=0.
    if (doplacn+time>yrmoda(2001,12,31)-yrmoda(1899,12,31)) reint=0.

    execute.

    * Create weights to adjust for time dependent censoring.

    * Double up the cases with no re-intervention and weight by probability of re-intervention.
    compute caseid=$casenum.
    execute.

    * Look up probability of reattendance after last visit.
    sort cases by interval.
    match files file=*
     /table='D:\Longitudinal Data\Rep_reattendance probability.sav'
     /by interval.
    compute recid=1.
    sort cases by caseid.
    save outfile='d:\temp1.sav'.

    select if (reint=0).
    compute recid=2.
    sort cases by caseid.
    save outfile='d:\temp2.sav'.

    add files file='d:\temp1.sav'
     /file='d:\temp2.sav'
     /by caseid .

    Compute weight = 1000.
    Do if (reint=0).
    if (recid=1) weight = 1000-trunc(1000*prreatt+0.5).
    if (recid=2) weight = trunc(1000*prreatt+0.5).
    end if.

    save outfile='d:\temp3.sav'.

    get file='d:\temp3.sav'.

    ***************************************************.
    * Compare and calculate survival function and statistics, and standard errors.

    * Assume censored cases are lost at one year (365 days) after last visit.

    * recalculate time if reint=0.
    * Impute effective time of censoring.
    Do if (reint=0).
    If (recid=1) time=min(lastvn+365,yrmoda(2001,12,31)-yrmoda(1899,12,31))- doplacn.
    if (recid=2) time=yrmoda(2001,12,31)-yrmoda(1899,12,31)- doplacn.
    end if.

    save outfile= 'd:\temp4.sav'.
    get file='d:\temp4.sav'.

    *  Test 1 Adjusted KM, but increased weights give incorrect confidence intervals.

    weight by weight.

    KM
      time  /STATUS=reint(1)
        /PRINT MEAN
        /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95)
        /Save survival (survadj).

# Kaplan-Meier

Notes
Output Created | 28-MAY-2003 12:18:13
Comments |
Input | Data | d:\temp4.sav
Filter | <none>
Weight | WEIGHT
Split File | <none>
N of Rows in Working Data File | 814068
Syntax | KM time /STATUS=reint(1) /PRINT MEAN /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95) /Save survival (survadj).
Resources | Elapsed Time | 0:01:59.25

     Survival Analysis for TIME

              Survival Time    Standard Error   95% Confidence Interval

     Mean:      2553.18                 .09     (  2553.01,   2553.36 )
     (Limited to   4017.0 )
     Median:    3061.00                 .46     (  3060.09,   3061.91 )

                                                                   Percentiles

                         45.00        50.00        55.00        60.00        65.00        70.00        75.00        80.00        85.00

     Value             3808.00      3061.00      2486.00      1997.00      1593.00      1255.00       967.00       723.00       515.00
     Standard Error        .            .46          .33          .25          .19          .15          .12          .09          .07

                                                                   Percentiles

                         90.00        95.00

     Value              332.00       175.00
     Standard Error        .05          .

    >Warning # 3211
    >On at least one case, the value of the weight variable was zero, negative,
    >or missing.  Such cases are invisible to statistical procedures and graphs
    >which need positively weighted cases, but remain on the file and are
    >processed by non-statistical facilities such as LIST and SAVE.

    * Test 2 Unadjusted KM, for confidence intervals.

    weight off.

    * Censor all non-returners.
    select if (recid=1).

    KM
      time  /STATUS=reint(1)
       /PRINT MEAN
       /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95)
       /Save survival (survunad).

# Kaplan-Meier

Notes
Output Created | 28-MAY-2003 12:20:18
Comments |
Input | Data | d:\temp4.sav
Filter | <none>
Weight | <none>
Split File | <none>
N of Rows in Working Data File | 485400
Syntax | KM time /STATUS=reint(1) /PRINT MEAN /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95) /Save survival (survunad).
Resources | Elapsed Time | 0:01:32.18

     Survival Analysis for TIME

              Survival Time    Standard Error   95% Confidence Interval

     Mean:      2498.92                2.89     (  2493.26,   2504.58 )
     (Limited to   4017.0 )
     Median:    2824.00               12.00     (  2800.48,   2847.52 )

                                                                   Percentiles

                         45.00        50.00        55.00        60.00        65.00        70.00        75.00        80.00        85.00

     Value             3416.00      2824.00      2320.00      1890.00      1519.00      1211.00       941.00       710.00       509.00
     Standard Error        .          12.00         9.20         7.24         5.66         4.44         3.50         2.70         2.04

                                                                   Percentiles

                         90.00        95.00

     Value              332.00       175.00
     Standard Error       1.47          .

    * Test 3 Censoring using randomisation instead of weighting.

    * Censor the proportion prreatt of all non-returners at the end of the observation period (instead of earlier).

    Set seed = 10102002.
    compute random=uniform(1).
    if (reint=0 and random<prreatt) time=yrmoda(2001,12,31)-yrmoda(1899,12,31)- doplacn.

    KM
      time  /STATUS=reint(1)
       /PRINT MEAN
       /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95)
       /Save survival (survrand).

# Kaplan-Meier

Notes
Output Created | 28-MAY-2003 12:21:50
Comments |
Input | Data | d:\temp4.sav
Filter | <none>
Weight | <none>
Split File | <none>
N of Rows in Working Data File | 485400
Syntax | KM time /STATUS=reint(1) /PRINT MEAN /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95) /Save survival (survrand).
Resources | Elapsed Time | 0:03:16.68

     Survival Analysis for TIME

              Survival Time    Standard Error   95% Confidence Interval

     Mean:      2552.64                2.82     (  2547.11,   2558.16 )
     (Limited to   4017.0 )
     Median:    3059.00               14.67     (  3030.25,   3087.75 )

                                                                   Percentiles

                         45.00        50.00        55.00        60.00        65.00        70.00        75.00        80.00        85.00

     Value             3807.00      3059.00      2485.00      1996.00      1592.00      1254.00       966.00       723.00       515.00
     Standard Error        .          14.67        10.42         7.99         6.11         4.74         3.67         2.80         2.10

                                                                   Percentiles

                         90.00        95.00

     Value              332.00       175.00
     Standard Error       1.49          .

    save outfile='D:\temp5.sav'.

    * Rerun to test robustness.

    get file='d:\temp4.sav'.
    select if (recid=1).
    compute random=uniform(1).
    if (reint=0 and random<prreatt) time=yrmoda(2001,12,31)-yrmoda(1899,12,31)- doplacn.

    KM
      time  /STATUS=reint(1)
       /PRINT MEAN
       /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95).

# Kaplan-Meier

Notes
Output Created | 28-MAY-2003 12:26:36
Comments |
Input | Data | d:\temp4.sav
Filter | <none>
Weight | <none>
Split File | <none>
N of Rows in Working Data File | 485400
Syntax | KM time /STATUS=reint(1) /PRINT MEAN /PERCENTILES = (45,50,55,60,65,70,75,80,85,90,95).
Resources | Elapsed Time | 0:00:43.02

     Survival Analysis for TIME

              Survival Time    Standard Error   95% Confidence Interval

     Mean:      2552.92                2.82     (  2547.40,   2558.44 )
     (Limited to   4017.0 )
     Median:    3059.00               14.67     (  3030.26,   3087.74 )

                                                                   Percentiles

                         45.00        50.00        55.00        60.00        65.00        70.00        75.00        80.00        85.00

     Value             3807.00      3059.00      2486.00      1997.00      1593.00      1255.00       967.00       723.00       515.00
     Standard Error        .          14.67        10.42         7.98         6.11         4.75         3.68         2.80         2.10

                                                                   Percentiles

                         90.00        95.00

     Value              332.00       175.00
     Standard Error       1.49          .

    **************************************************************************************.

    * Create an Excel file for comparing survivor functions.
    get file='d:\temp5.sav'.

    weight off.

    select if (reint=1).
    sort cases by time.
    aggregate outfile=*
     /presorted
     /break time
     / survadj survunad survrand
    = first( survadj survunad survrand ).

    SAVE TRANSLATE OUTFILE='I:\Research Projects\longevity\Rep_1991 to 2001 adj unadj and rand.xls'
      /TYPE=XLS /MAP /REPLACE /FIELDNAMES.

    Data written to I:\Research Projects\longevity\Rep_1991 to 2001 adj unadj and rand.xls.
    4 variables and 3775 cases written to range: SPSS.
    Variable: TIME       Type: Number   Width:  8   Dec: 2
    Variable: SURVADJ    Type: Number   Width: 10   Dec: 5
    Variable: SURVUNAD   Type: Number   Width: 10   Dec: 5
    Variable: SURVRAND   Type: Number   Width: 10   Dec: 5

    * Create an Excel file for calculation of standard errors.

    get file='d:\temp4.sav'.

    weight by weight.

    sort cases by time reint.

    save outfile = 'd:\temp6.sav'.
    get file='d:\temp6.sav'.

    compute censored=0.
    if (reint=0) censored=1.

    aggregate outfile=*
     /presorted
     /break time
     /reints censored =sum(reint censored ).

    >Warning # 3211
    >On at least one case, the value of the weight variable was zero, negative,
    >or missing.  Such cases are invisible to statistical procedures and graphs
    >which need positively weighted cases, but remain on the file and are
    >processed by non-statistical facilities such as LIST and SAVE.

    compute reints=reints/1000.
    compute censored = censored/1000.

    SAVE TRANSLATE OUTFILE='I:\Research Projects\longevity\Rep_1991 to 2001 times for se.xls'
      /TYPE=XLS /MAP /REPLACE /FIELDNAMES.

    Data written to I:\Research Projects\longevity\Rep_1991 to 2001 times for se.xls.
    3 variables and 4010 cases written to range: SPSS.
    Variable: TIME       Type: Number   Width:  8   Dec: 2
    Variable: REINTS     Type: Number   Width:  8   Dec: 2
    Variable: CENSORED   Type: Number   Width:  8   Dec: 2