each covariate and the intervention selection is pre-specified and in a fixed form. However, one main challenge of the a priori modeling approach is that it might not provide the optimal balance between treatment and control groups.
To build an a priori model for propensity score estimation in SAS, we can use either PROC PSMATCH or PROC LOGISTIC as shown in Program 4.5. In both cases, the input data set is a one observation per patient data set containing the treatment and baseline covariates (from the simulated REFLECTIONS study – see Chapter 3). Also, in both cases the code will produce an output data set containing the original data set with the additional estimated propensity score for each patient (_ps_).
Program 4.5: Propensity Score Estimation: A priori Logistic Regression
PROC PSMATCH DATA=REFL2 REGION=ALLOBS;
CLASS COHORT GENDER RACE DR_RHEUM DR_PRIMCARE;
PSMODEL COHORT(TREATED=’OPIOID’)= GENDER RACE AGE BMI_B BPIINTERF_B BPIPAIN_B
CPFQ_B FIQ_B GAD7_B ISIX_B PHQ8_B PHYSICALSYMP_B SDS_B DR_RHEUM
DR_PRIMCARE;
OUTPUT OUT=PS PS=_PS_;
RUN;
PROC LOGISTIC DATA=REFL2;
CLASS COHORT GENDER RACE DR_RHEUM DR_PRIMCARE;
MODEL COHORT = GENDER RACE AGE BMI_B BPIINTERF_B BPIPAIN_B CPFQ_B FIQ_B GAD7_B
ISIX_B PHQ8_B PHYSICALSYMP_B SDS_B DR_RHEUM DR_PRIMCARE;
OUTPUT OUT=PS PREDICTED=PS;
RUN;
Before building a logistic model in SAS, we suggest examining the distribution of the intervention indicator at each level of the categorical variable to rule out the possibility of “complete separation” (or “perfect prediction”), which means that for subjects at some level of some categorical variable, they would all receive one intervention but not the other. Complete separation can occur for several reasons and one common example is when using several categorical variables whose categories are coded by indicators. When the logistic regression model is fit, the estimate of the regression coefficients
WARNING: There is a complete separation of data points. The maximum likelihood estimate does not exist.
WARNING: The LOGISTIC procedure continues in spite of the above warning. Results shown are based on the last maximum likelihood iteration. Validity of the model fit is questionable.
Notice that SAS will continue to finish the computation despite issuing warning messages. However, the estimate of such s are incorrect, and so are the estimated propensity scores. If after examining the intervention distribution at each level of the categorical variables complete separation is found, then efforts should be made to address this issue. One possible solution is to collapse the categorical variable causing the problem. That is, combine the different outcome categories such that the complete separation no longer exists. Another possible solution is to use Firth logistic regression. It uses a penalized likelihood estimation method. Firth bias-correction is considered an ideal solution to the separation issue for logistic regression (Heinze and Schemper, 2002). In PROC LOGISTIC, we can add an option to run the Firth logistic regression as shown in Program 4.6.
Program 4.6: Firth Logistic Regression
PROC LOGISTIC DATA=REFL2;
CLASS COHORT GENDER RACE DR_RHEUM DR_PRIMCARE;
MODEL COHORT = GENDER RACE DR_RHEUM DR_PRIMCARE BPIInterf_B BPIPain_B
CPFQ_B FIQ_B GAD7_B ISIX_B PHQ8_B PhysicalSymp_B SDS_B / FIRTH;
OUTPUT OUT=PS PREDICTED=PS;
RUN;
Automatic Parametric Model Selection
A second parametric approach to estimate the propensity score is the use of an automated model building process that ensures the balance of confounders between interventions. The idea originated from Rubin and Rosenbaum (Rubin and Rosenbaum, 1984) and was later developed and proposed in detail by Dehejia and Wahba (Dehejia and Wahba, 1998, 1999). This approach also has broad application in other areas such as psychology (Stuart, 2003) and economics (Marco and Kopeinig, 2008). It is an iterative approach using logistic regression and we suggest using the following steps for implementation:
1. Estimate propensity scores using a logistic regression model that included the treatment indicator as the dependent variable and other measured covariates as explanatory variables. No interaction terms or high-order terms of those covariates are included at this step.
2. Order the estimated propensity score in step 1 from low to high, then divide the propensity scores into strata such that each stratum holds an approximately equal number of treated individuals. Some studies (Stuart, 2003) suggested five strata would be a reasonable choice to avoid too few comparison subjects within a stratum.
3. Evaluate the balance of the measured covariates and all of their two-way interaction-terms within each stratum from step 2. Balance can be quantified using the standardized bias (or standardized mean difference, see Chapter 5). For continuous covariates, the standardized bias is the difference in means of the covariate between the treated group and the comparison group divided by the standard deviation of the treatment group. For categorical covariates, at each level of the categorical covariate, the standardized bias is the difference in proportions of each level of the measured covariate divided by the standard deviation in the treatment group. To be more precise:
◦ For interactions between continuous covariates, create a new variable that is the product of the two variables.
◦ For interactions between a categorical and a continuous covariate, calculate the standardized difference per level of the categorical variable.
◦ For interactions between two categorical covariates A and B we calculate the following:
- For each level of A, the difference in proportions for each level of B divided by the standard deviation in the treatment group.
- For each level of B, the difference in proportions for each level of A divided by the standard deviation in the treatment group.
A covariate is considered “balanced” if the standardized bias is less than 0.25. Although some suggest a stricter threshold, < 0.1, which is preferred in propensity score-based analysis. If the standardized bias for an interaction is above 0.25 across two or more strata, then that interaction term is considered “imbalanced.”
4. After adding each interaction term to the model separately, keep the one that reduces the number of imbalanced interaction terms the most (in other words, improve the balance the most). Then fit all remaining interaction terms again separately, repeat steps 2 and 3, and again add the one that improves balance most. Repeat until no more improvement.
Program 4.9 provides the macros to implement the automatic model selection for estimating propensity score in SAS. For preparation, Program 4.7 provides a macro to automatically create binary indicator for all categorical variables specified in Programs 4.9 and 4.10. Notice that Program 4.7 will create (n-1) dummy binary variables for a categorical variable with n categories if this