3 Interaction without adjustment
fit_unadjusted <- mfpi(
Surv(survtime, cens) ~ rx + fp(age),
data = advanced_prostate_cancer,
family = "cox",
group_var = "rx",
interaction_vars = "age",
interaction_forms = c(age = "fp2"),
flex = "flex1",
ties = "breslow"
)
#>
#> ======================================================================
#> STEP 1: Selected Adjustment Model (MFP)
#> ======================================================================
#>
#> Model: Cox proportional hazards
#> Criterion: p-value
#>
#> df_initial select alpha searched_fp df_final power1
#> age 4 0.050 0.050 yes 2 3
#>
#>
#> Selected adjustment variables (1): age
#>
#> ======================================================================
#> STEP 2: Evaluating Interactions (flex = flex1)
#> ======================================================================
#>
#> Group variable: rx
#> Selection rule: p-value < 0.050
#> P-value adjustment: none
#>
#> ----------------------------------------------------------------------
#> [1/1] Variable: 'age' x 'rx'
#> Adjusting for: (none)
#> ----------------------------------------------------------------------
#> Type deviance_int deviance_diff df pvalue main_power int_power
#> FP2 3766.216 9.746 2 0.008 (3, 3) (3, 3), (3, 3)
#>
#>
#> ======================================================================
#> STEP 3: Interaction Summary
#> ======================================================================
#> Variable Type p-value selected
#> ----------------------------------------------------------------------
#> age FP2 0.008 Yes
#> ----------------------------------------------------------------------
#>
#> Done: 1 of 1 interactions selected.
print(fit_unadjusted)
#> =========================================
#> MFPI | n = 475 (events = 338) | flex1
#> =========================================
#>
#> Call:
#> mfpi(formula = Surv(survtime, cens) ~ rx + fp(age), data = advanced_prostate_cancer,
#> group_var = "rx", interaction_vars = "age", interaction_forms = c(age = "fp2"),
#> flex = "flex1", family = "cox", ties = "breslow")
#>
#> Model: Cox proportional hazards
#>
#>
#> Step 1 - Selected Adjustment Model (MFP):
#> -----------------------------------------
#> df_initial select alpha searched_fp df_final power1
#> age 4 0.050 0.050 yes 2 3
#>
#> Selected adjustment variables (1): age
#>
#> Step 2 - All Interaction Candidates:
#> ------------------------------------
#>
#> Grouping variable : rx
#> Selection criterion : p-value
#> Selection rule : p-value < 0.050
#> P-value adjustment : none
#>
#> variable type deviance_int deviance_diff df_int p-value fp_powers_main
#> age FP2 3766.216 9.746 2 0.008 (3, 3)
#>
#> FP powers by group level:
#> variable type Group Group level fp_powers_int
#> age FP2 rx Placebo (3, 3)
#> age FP2 rx Treatment (3, 3)
#>
#> Step 3 - Interaction Summary:
#> -----------------------------
#> variable type p-value selected
#> age FP2 0.008 YesWhat the arguments do:
fp(age)puts age in the Stage 1 adjustment model. It does not set the form of the interaction.group_varnames the grouping variable.interaction_varslists the continuous variables tested for interaction withgroup_var.interaction_forms = c(age = "fp2")models the interaction with an FP2 function. All 36 FP2 power pairs are compared. Variables not listed use FP1.flex = "flex1"selects the FP powers from the model without interaction and uses them in both groups (paper, Section 2.5).
How to read the output (paper, Section 4.2):
- Step 1 is the adjustment model. Age is selected as FP1 with power 3. Here age is the only adjustment variable, and it is removed from the adjustment set when its own interaction is tested.
- Step 2 is the interaction test. The FP2 powers (3, 3) are used in both groups. Adding the interaction reduces the deviance by 9.746 on 2 df.
- Step 3 summarises the decision: the interaction is selected (p = 0.008).
3.1 Plotting the interaction
plot() returns a list with one element per variable.
Each element contains one plot per comparison with the reference
group.
Unadjusted analysis: fitted age functions by group (left) and the Treatment-minus-Placebo log hazard ratio with 95% CI (right); paper, Figure 2.
The right panel shows the log hazard ratio for Treatment vs Placebo as a function of age, the difference function \(\Delta\) (paper, Section 2.3). The dashed line is the constant log hazard ratio from the model without interaction. Here, the benefit of treatment is largest in younger patients and decreases with age. The confidence intervals do not account for the selection of FP powers.
3.2 Predictions at chosen ages
type = "both" returns the fitted function for each group
and the difference between groups, all on the linear predictor
scale.
age_pred <- predict(
fit_unadjusted,
newdata = data.frame(age = c(50, 60, 70)),
terms = "age",
level = 0.95,
type = "both"
)
age_pred
#> =========================================================
#> MFPI prediction | term: age | scale: linear predictor
#> =========================================================
#>
#> Partial fitted functions:
#> age group estimate SE lower upper
#> 50 Placebo 0.460 0.322 -0.170 1.090
#> 60 Placebo -0.018 0.126 -0.266 0.230
#> 70 Placebo -0.154 0.080 -0.310 0.003
#> 50 Treatment -1.014 0.404 -1.805 -0.223
#> 60 Treatment -0.821 0.199 -1.212 -0.430
#> 70 Treatment -0.376 0.123 -0.618 -0.134
#>
#> Differences:
#> age contrast group reference estimate SE lower upper
#> 50 Treatment-Placebo Treatment Placebo -1.474 0.517 -2.488 -0.460
#> 60 Treatment-Placebo Treatment Placebo -0.802 0.235 -1.263 -0.342
#> 70 Treatment-Placebo Treatment Placebo -0.223 0.144 -0.505 0.060
#>
#> Note: Partial fitted functions comprise the intercept (when present), group main effect,
#> and group-specific FP function estimated from the adjusted interaction model.
#> Differences are comparison-group minus reference-group partial fitted functions.The differences are age-specific log hazard ratios. Exponentiating them gives hazard ratios of about 0.23, 0.45 and 0.80 at ages 50, 60 and 70 (paper, Section 4.2.2).
