Three Age-Group Analysis

Library and Data

library(survival)
library(tidyverse)
library(penalized)
library(kableExtra)
library(stringr)
library(patchwork)

source("../R/utils.R")

library(readxl)
raw_data <- read_excel("../data/Aim3-early-adversity-subjects-F-survival-2021-01-07.xlsx")

Raw data here refers to the larger dataset.

EDA

Below is the Kaplan-Meier estimation of the survival function on the raw data.

km <- survfit(Surv(raw_data$age, raw_data$dead) ~ 1, conf.type = "log-log")
ggplot() +
  geom_line(aes(x = km$time, y = km$surv)) +
  geom_vline(aes(xintercept = 1), color = "red", linetype = "dashed") +
  geom_vline(aes(xintercept = 5), color = "red", linetype = "dashed") +
  labs(x = expression("Time (year)"), y = expression("Survival"),
       title = "Kaplan-Meier estimate") +
  scale_x_continuous(breaks = seq(0, 25, 1)) +
  theme_classic()

Analysis Process

We conduct four models on three datasets derived from the raw data.

Four models:

Cox model of survival data versus all the main effects (six adversities)
Cox model of survival data versus main effects and interaction terms selected by lasso
Cox model of survival data versus main effects and interaction terms selected by stepwise BIC
Cox model of survival data versus adv_mom, sum of other adversities adv_other_cumulative and their interaction term

Three datasets:

age_below_1: Include all the samples. For those whose age > 1, use age = 1, and dead = FALSE
age_above_1_below_5: Include samples of age > 1. For those whose age > 5, use age = 5, and dead = FALSE
age_above_5: Include samples of age > 5.

Note that for the age_below_1 dataset, we exclude adv_sib from all models.

age_below_1 <- raw_data %>%
  select(starts_with("adv_"), age, dead, -adv_cumulative) %>%
  mutate(dead = if_else(age > 1, FALSE, dead)) %>%
  mutate(age = if_else(age > 1, 1, age))

age_above_1_below_5 <- raw_data %>%
  filter(age >= 1) %>%
  select(starts_with("adv_"), age, dead, -adv_cumulative) %>%
  mutate(dead = if_else(age > 5, FALSE, dead)) %>%
  mutate(age = if_else(age > 5, 5, age))

age_above_5 <- raw_data %>%
  filter(age >= 5) %>%
  select(starts_with("adv_"), age, dead, -adv_cumulative)

Below is the summary of three datasets. The number of events refers to the number of uncensored death.

data.frame(
  nobs = c(nrow(age_below_1), nrow(age_above_1_below_5), nrow(age_above_5)),
  nevt = c(sum(age_below_1$dead), sum(age_above_1_below_5$dead), sum(age_above_5$dead))
) %>%
  `row.names<-`(c("Age 1-", "Age 1+, 5-", "Age 5+")) %>%
  kable(row.names = TRUE,
        col.names = c("Number of Samples", "Number of Events")) %>%
  kable_styling(full_width = FALSE, position = "center", font_size = 18)

	Number of Samples	Number of Events
Age 1-	562	120
Age 1+, 5-	440	92
Age 5+	277	110

Cox Model - Main Effects

run_full_model <- function(data, below_one = FALSE) {
  varnames = names(data)[startsWith(names(data), "adv_")]
  if (below_one) {
    varnames = varnames[varnames != "adv_sib"]
  }
  full_formula = paste("Surv(age, dead) ~", paste0(varnames, collapse = "+ "))
  full_model = coxph(
    as.formula(full_formula),
    data = data
  )
  return (full_model)
}

Age below 1

full_age_below_1 <- run_full_model(age_below_1, below_one = TRUE)
summary(full_age_below_1)

## Call:
## coxph(formula = as.formula(full_formula), data = data)
## 
##   n= 562, number of events= 120 
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)  
## adv_densityTRUE  -0.4879    0.6139   0.2529 -1.929   0.0537 .
## adv_rainTRUE     -0.2429    0.7844   0.3179 -0.764   0.4449  
## adv_momTRUE       0.4337    1.5429   0.1946  2.229   0.0258 *
## adv_mom_rankTRUE  0.1504    1.1624   0.2084  0.722   0.4704  
## adv_mom_dsiTRUE  -0.1680    0.8453   0.2214 -0.759   0.4480  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                  exp(coef) exp(-coef) lower .95 upper .95
## adv_densityTRUE     0.6139     1.6289    0.3740     1.008
## adv_rainTRUE        0.7844     1.2749    0.4206     1.463
## adv_momTRUE         1.5429     0.6481    1.0538     2.259
## adv_mom_rankTRUE    1.1624     0.8603    0.7726     1.749
## adv_mom_dsiTRUE     0.8453     1.1829    0.5477     1.305
## 
## Concordance= 0.575  (se = 0.025 )
## Likelihood ratio test= 11.53  on 5 df,   p=0.04
## Wald test            = 11.23  on 5 df,   p=0.05
## Score (logrank) test = 11.43  on 5 df,   p=0.04

Age above 1 and below 5

full_age_above_1_below_5 <- run_full_model(age_above_1_below_5)
summary(full_age_above_1_below_5)

## Call:
## coxph(formula = as.formula(full_formula), data = data)
## 
##   n= 440, number of events= 92 
## 
##                      coef exp(coef) se(coef)      z Pr(>|z|)    
## adv_densityTRUE   0.04178   1.04266  0.24461  0.171 0.864388    
## adv_rainTRUE     -0.36927   0.69124  0.33866 -1.090 0.275537    
## adv_momTRUE       0.78781   2.19858  0.21648  3.639 0.000274 ***
## adv_sibTRUE      -0.08060   0.92256  0.27909 -0.289 0.772730    
## adv_mom_rankTRUE  0.22330   1.25020  0.24313  0.918 0.358394    
## adv_mom_dsiTRUE  -0.61402   0.54117  0.27299 -2.249 0.024497 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                  exp(coef) exp(-coef) lower .95 upper .95
## adv_densityTRUE     1.0427     0.9591    0.6455    1.6841
## adv_rainTRUE        0.6912     1.4467    0.3559    1.3424
## adv_momTRUE         2.1986     0.4548    1.4384    3.3606
## adv_sibTRUE         0.9226     1.0839    0.5339    1.5943
## adv_mom_rankTRUE    1.2502     0.7999    0.7763    2.0134
## adv_mom_dsiTRUE     0.5412     1.8478    0.3169    0.9241
## 
## Concordance= 0.619  (se = 0.029 )
## Likelihood ratio test= 18.38  on 6 df,   p=0.005
## Wald test            = 18.7  on 6 df,   p=0.005
## Score (logrank) test = 19.48  on 6 df,   p=0.003

Age above 5

full_age_above_5 <- run_full_model(age_above_5)
summary(full_age_above_5)

## Call:
## coxph(formula = as.formula(full_formula), data = data)
## 
##   n= 277, number of events= 110 
## 
##                     coef exp(coef) se(coef)      z Pr(>|z|)   
## adv_densityTRUE   0.1406    1.1510   0.2672  0.526  0.59864   
## adv_rainTRUE     -0.2902    0.7481   0.3099 -0.936  0.34903   
## adv_momTRUE       0.7398    2.0956   0.2257  3.278  0.00105 **
## adv_sibTRUE       0.2405    1.2719   0.2469  0.974  0.33003   
## adv_mom_rankTRUE  0.2448    1.2773   0.2402  1.019  0.30822   
## adv_mom_dsiTRUE  -0.2983    0.7421   0.2320 -1.286  0.19852   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                  exp(coef) exp(-coef) lower .95 upper .95
## adv_densityTRUE     1.1510     0.8688    0.6818     1.943
## adv_rainTRUE        0.7481     1.3367    0.4075     1.373
## adv_momTRUE         2.0956     0.4772    1.3464     3.262
## adv_sibTRUE         1.2719     0.7862    0.7839     2.064
## adv_mom_rankTRUE    1.2773     0.7829    0.7977     2.045
## adv_mom_dsiTRUE     0.7421     1.3476    0.4709     1.169
## 
## Concordance= 0.616  (se = 0.035 )
## Likelihood ratio test= 13.93  on 6 df,   p=0.03
## Wald test            = 14.83  on 6 df,   p=0.02
## Score (logrank) test = 15.31  on 6 df,   p=0.02

Cox Model - Lasso

get_penalized_covariates <- function(data, below_one=FALSE) {
  data = select(data, starts_with("adv_"))
  if (below_one) {
    data = select(data, -adv_sib)
  }
  penalized_covariates <- model.matrix(
    as.formula("~ . + .:."),
    data = data.frame(data)
  )
  penalized_covariates <- penalized_covariates[, colnames(penalized_covariates) != "(Intercept)"]
  return (penalized_covariates)
}

get_opt_lambda_of_lasso <- function(data, penalized_covariates, seed = 42, fold = 5) {
  set.seed(seed)
  opt <- optL1(
    Surv(data$age, data$dead),
    penalized = penalized_covariates,
    standardize = T,
    fold = fold,
    trace = F
  )
  return (opt$lambda)
}

get_nonzero_varnames_of_opt_lasso <- function(data, penalized_covariates, opt_lambda) {
  final = penalized(
    Surv(data$age, data$dead),
    penalized = penalized_covariates,
    standardize = T,
    lambda1 = opt_lambda,
    trace = F
  )
  return (names(final@penalized)[final@penalized != 0])
}

get_varnames_in_final_model <- function(nonzero_varnames) {
  # append main effects select by interaction terms
  ans = c()
  for (varname in nonzero_varnames) {
    if (str_detect(varname, ":")) {
      ans = c(ans, str_split(varname, ":")[[1]])
    }
    ans = c(ans, varname)
  }
  return (unique(ans))
}

run_cox <- function(data, varnames, penalized_covariates) {
  formula = paste("Surv(age, dead) ~", paste0(varnames, collapse = "+ "))
  cox <- coxph(
    as.formula(formula),
    data = data.frame(
      cbind(penalized_covariates,
            select(data, c("age", "dead"))
      ))
  )
  return (cox)
}

Age below 5

penalized_below_1 <- get_penalized_covariates(age_below_1, below_one = TRUE)
opt_lambda_below_1 <- get_opt_lambda_of_lasso(
  data = age_below_1,
  penalized_covariates = penalized_below_1
)

nonzero_vars_below_1 <- get_nonzero_varnames_of_opt_lasso(
  data = age_below_1,
  penalized_covariates = penalized_below_1,
  opt_lambda = opt_lambda_below_1
)

vars_below_1 <- get_varnames_in_final_model(nonzero_vars_below_1)
cox_below_1 <- run_cox(
  data = age_below_1,
  varnames = vars_below_1,
  penalized_covariates = penalized_below_1
)

summary(cox_below_1)

## Call:
## coxph(formula = as.formula(formula), data = data.frame(cbind(penalized_covariates, 
##     select(data, c("age", "dead")))))
## 
##   n= 562, number of events= 120 
## 
##                    coef exp(coef) se(coef)      z Pr(>|z|)  
## adv_densityTRUE -0.5308    0.5881   0.2501 -2.123   0.0338 *
## adv_momTRUE      0.4401    1.5528   0.1937  2.272   0.0231 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                 exp(coef) exp(-coef) lower .95 upper .95
## adv_densityTRUE    0.5881      1.700    0.3602    0.9601
## adv_momTRUE        1.5528      0.644    1.0623    2.2698
## 
## Concordance= 0.565  (se = 0.022 )
## Likelihood ratio test= 9.94  on 2 df,   p=0.007
## Wald test            = 9.66  on 2 df,   p=0.008
## Score (logrank) test = 9.85  on 2 df,   p=0.007

Age above 1 and below 5

penalized_above_1_below_5 <- get_penalized_covariates(age_above_1_below_5)
opt_lambda_above_1_below_5 <- get_opt_lambda_of_lasso(
  data = age_above_1_below_5,
  penalized_covariates = penalized_above_1_below_5
)
nonzero_vars_above_1_below_5 <- get_nonzero_varnames_of_opt_lasso(
  data = age_above_1_below_5,
  penalized_covariates = penalized_above_1_below_5,
  opt_lambda = opt_lambda_above_1_below_5
)

vars_above_1_below_5 <- get_varnames_in_final_model(nonzero_vars_above_1_below_5)
cox_above_1_below_5 <- run_cox(
  data = age_above_1_below_5,
  varnames = vars_above_1_below_5,
  penalized_covariates = penalized_above_1_below_5
)

summary(cox_above_1_below_5)

## Call:
## coxph(formula = as.formula(formula), data = data.frame(cbind(penalized_covariates, 
##     select(data, c("age", "dead")))))
## 
##   n= 440, number of events= 92 
## 
##                                 coef exp(coef) se(coef)      z Pr(>|z|)  
## adv_momTRUE                   0.5600    1.7506   0.2585  2.166   0.0303 *
## adv_mom_dsiTRUE              -0.6434    0.5255   0.2743 -2.346   0.0190 *
## adv_mom_rankTRUE             -0.1233    0.8840   0.3350 -0.368   0.7128  
## adv_momTRUE:adv_mom_rankTRUE  0.8102    2.2484   0.4958  1.634   0.1022  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                              exp(coef) exp(-coef) lower .95 upper .95
## adv_momTRUE                     1.7506     0.5712    1.0547    2.9058
## adv_mom_dsiTRUE                 0.5255     1.9030    0.3070    0.8996
## adv_mom_rankTRUE                0.8840     1.1313    0.4585    1.7044
## adv_momTRUE:adv_mom_rankTRUE    2.2484     0.4448    0.8509    5.9412
## 
## Concordance= 0.612  (se = 0.029 )
## Likelihood ratio test= 19.7  on 4 df,   p=6e-04
## Wald test            = 21.23  on 4 df,   p=3e-04
## Score (logrank) test = 22.39  on 4 df,   p=2e-04

Age above 5

penalized_above_5 <- get_penalized_covariates(age_above_5)
opt_lambda_above_5 <- get_opt_lambda_of_lasso(
  data = age_above_5,
  penalized_covariates = penalized_above_5
)
nonzero_vars_above_5 <- get_nonzero_varnames_of_opt_lasso(
  data = age_above_5,
  penalized_covariates = penalized_above_5,
  opt_lambda = opt_lambda_above_5
)

vars_above_5 <- get_varnames_in_final_model(nonzero_vars_above_5)
cox_above_5 <- run_cox(
  data = age_above_5,
  varnames = vars_above_5,
  penalized_covariates = penalized_above_5
)

summary(cox_above_5)

## Call:
## coxph(formula = as.formula(formula), data = data.frame(cbind(penalized_covariates, 
##     select(data, c("age", "dead")))))
## 
##   n= 277, number of events= 110 
## 
##                             coef exp(coef) se(coef)      z Pr(>|z|)   
## adv_momTRUE               0.8824    2.4166   0.2720  3.244  0.00118 **
## adv_rainTRUE              0.3574    1.4296   0.3335  1.072  0.28384   
## adv_sibTRUE               0.2025    1.2245   0.2935  0.690  0.49020   
## adv_momTRUE:adv_rainTRUE -1.8708    0.1540   0.8200 -2.282  0.02251 * 
## adv_momTRUE:adv_sibTRUE   0.3785    1.4601   0.5123  0.739  0.45998   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                          exp(coef) exp(-coef) lower .95 upper .95
## adv_momTRUE                  2.417     0.4138   1.41804    4.1184
## adv_rainTRUE                 1.430     0.6995   0.74362    2.7486
## adv_sibTRUE                  1.224     0.8167   0.68885    2.1765
## adv_momTRUE:adv_rainTRUE     0.154     6.4938   0.03087    0.7682
## adv_momTRUE:adv_sibTRUE      1.460     0.6849   0.53495    3.9855
## 
## Concordance= 0.623  (se = 0.031 )
## Likelihood ratio test= 18.96  on 5 df,   p=0.002
## Wald test            = 22.11  on 5 df,   p=5e-04
## Score (logrank) test = 24.41  on 5 df,   p=2e-04

Cox Model - Stepwise BIC

Age below 1

full_bic_below_1 <- coxph(
  as.formula(Surv(age, dead) ~ . + .:.),
  data = select(age_below_1, -adv_sib)
)

step_bic_below_1 <- step(
  full_bic_below_1,
  k = log(nrow(age_below_1)),
  #scope = list(lower= ~ adv_density + adv_rain + adv_mom + adv_mom_rank + 
  #  adv_mom_dsi),
  trace = F
)

summary(step_bic_below_1)

## Call:  coxph(formula = Surv(age, dead) ~ 1, data = select(age_below_1, 
##     -adv_sib))
## 
## Null model
##   log likelihood= -746.0115 
##   n= 562

Age above 1, below 5

full_bic_above_1_below_5 <- coxph(
  as.formula(Surv(age, dead) ~ . + .:.),
  data = age_above_1_below_5
)

step_bic_above_1_below_5 <- step(
  full_bic_above_1_below_5,
  k = log(nrow(age_above_1_below_5)),
  #scope = list(lower= ~ adv_density + adv_rain + adv_mom + adv_mom_rank + 
  #  adv_mom_dsi + adv_sib),
  trace = F
)

summary(step_bic_above_1_below_5)

## Call:
## coxph(formula = Surv(age, dead) ~ adv_mom, data = age_above_1_below_5)
## 
##   n= 440, number of events= 92 
## 
##               coef exp(coef) se(coef)     z Pr(>|z|)    
## adv_momTRUE 0.7654    2.1498   0.2160 3.543 0.000395 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##             exp(coef) exp(-coef) lower .95 upper .95
## adv_momTRUE      2.15     0.4652     1.408     3.283
## 
## Concordance= 0.579  (se = 0.025 )
## Likelihood ratio test= 11.52  on 1 df,   p=7e-04
## Wald test            = 12.55  on 1 df,   p=4e-04
## Score (logrank) test = 13.18  on 1 df,   p=3e-04

Age above 5

full_bic_above_5 <- coxph(
  as.formula(Surv(age, dead) ~ . + .:.),
  data = age_above_5
)

step_bic_above_5 <- step(
  full_bic_above_5,
  k = log(nrow(age_above_5)),
  #scope = list(lower= ~ adv_density + adv_rain + adv_mom + adv_mom_rank + 
  #  adv_mom_dsi + adv_sib),
  trace = F
)

summary(step_bic_above_5)

## Call:
## coxph(formula = Surv(age, dead) ~ adv_rain + adv_mom + adv_rain:adv_mom, 
##     data = age_above_5)
## 
##   n= 277, number of events= 110 
## 
##                             coef exp(coef) se(coef)      z Pr(>|z|)    
## adv_rainTRUE              0.3199    1.3769   0.3288  0.973    0.331    
## adv_momTRUE               0.9748    2.6507   0.2342  4.162 3.16e-05 ***
## adv_rainTRUE:adv_momTRUE -1.9110    0.1479   0.8144 -2.346    0.019 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                          exp(coef) exp(-coef) lower .95 upper .95
## adv_rainTRUE                1.3769     0.7263   0.72289     2.623
## adv_momTRUE                 2.6507     0.3773   1.67490     4.195
## adv_rainTRUE:adv_momTRUE    0.1479     6.7596   0.02998     0.730
## 
## Concordance= 0.609  (se = 0.027 )
## Likelihood ratio test= 16.76  on 3 df,   p=8e-04
## Wald test            = 18.74  on 3 df,   p=3e-04
## Score (logrank) test = 20.24  on 3 df,   p=2e-04

Cox model - Cumulative

run_cox_with_cumulative <- function(data) {
  data = data %>%
    mutate(adv_other_cumulative = adv_sib+adv_rain+adv_mom_dsi+adv_mom_rank+adv_density) %>%
    select(age, dead, adv_mom, adv_other_cumulative)
  cox = coxph(
    formula = as.formula("Surv(age, dead) ~ . + adv_mom:adv_other_cumulative"),
    data = data
  )
  return (cox)
}

Age below 1

cox_below_1 <- run_cox_with_cumulative(age_below_1)
summary(cox_below_1)

## Call:
## coxph(formula = as.formula("Surv(age, dead) ~ . + adv_mom:adv_other_cumulative"), 
##     data = data)
## 
##   n= 562, number of events= 120 
## 
##                                     coef exp(coef) se(coef)      z Pr(>|z|)    
## adv_momTRUE                       0.8931    2.4427   0.3214  2.779  0.00545 ** 
## adv_other_cumulative              0.4640    1.5904   0.1113  4.168 3.07e-05 ***
## adv_momTRUE:adv_other_cumulative -0.3634    0.6953   0.1971 -1.843  0.06530 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                                  exp(coef) exp(-coef) lower .95 upper .95
## adv_momTRUE                         2.4427     0.4094    1.3011     4.586
## adv_other_cumulative                1.5904     0.6288    1.2786     1.978
## adv_momTRUE:adv_other_cumulative    0.6953     1.4381    0.4725     1.023
## 
## Concordance= 0.618  (se = 0.025 )
## Likelihood ratio test= 21.39  on 3 df,   p=9e-05
## Wald test            = 20.9  on 3 df,   p=1e-04
## Score (logrank) test = 21.32  on 3 df,   p=9e-05

Age above 1 and below 5

cox_above_1_below_5 <- run_cox_with_cumulative(age_above_1_below_5)
summary(cox_above_1_below_5)

## Call:
## coxph(formula = as.formula("Surv(age, dead) ~ . + adv_mom:adv_other_cumulative"), 
##     data = data)
## 
##   n= 440, number of events= 92 
## 
##                                      coef exp(coef) se(coef)      z Pr(>|z|)  
## adv_momTRUE                       0.70391   2.02164  0.32593  2.160   0.0308 *
## adv_other_cumulative             -0.17439   0.83997  0.16013 -1.089   0.2761  
## adv_momTRUE:adv_other_cumulative  0.07605   1.07902  0.24430  0.311   0.7556  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                                  exp(coef) exp(-coef) lower .95 upper .95
## adv_momTRUE                          2.022     0.4946    1.0673     3.829
## adv_other_cumulative                 0.840     1.1905    0.6137     1.150
## adv_momTRUE:adv_other_cumulative     1.079     0.9268    0.6685     1.742
## 
## Concordance= 0.602  (se = 0.031 )
## Likelihood ratio test= 13.04  on 3 df,   p=0.005
## Wald test            = 13.82  on 3 df,   p=0.003
## Score (logrank) test = 14.62  on 3 df,   p=0.002

Age above 5

cox_above_5 <- run_cox_with_cumulative(age_above_5)
summary(cox_above_5)

## Call:
## coxph(formula = as.formula("Surv(age, dead) ~ . + adv_mom:adv_other_cumulative"), 
##     data = data)
## 
##   n= 277, number of events= 110 
## 
##                                      coef exp(coef) se(coef)      z Pr(>|z|)   
## adv_momTRUE                       0.98709   2.68342  0.35927  2.748    0.006 **
## adv_other_cumulative              0.09773   1.10267  0.14777  0.661    0.508   
## adv_momTRUE:adv_other_cumulative -0.29894   0.74161  0.29290 -1.021    0.307   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
##                                  exp(coef) exp(-coef) lower .95 upper .95
## adv_momTRUE                         2.6834     0.3727    1.3270     5.426
## adv_other_cumulative                1.1027     0.9069    0.8254     1.473
## adv_momTRUE:adv_other_cumulative    0.7416     1.3484    0.4177     1.317
## 
## Concordance= 0.592  (se = 0.031 )
## Likelihood ratio test= 9.87  on 3 df,   p=0.02
## Wald test            = 11.14  on 3 df,   p=0.01
## Score (logrank) test = 11.67  on 3 df,   p=0.009

Discussion

For the full model with main effects of all six adversities, group age_below_1 and age_above_5 behave similarly in terms of the significant variable adv_mom. Group age_above_1_below_5 has one more significant relationship with adv_dsi, and the significance level of adv_mom is much lower, indicating the hazard might work differently within this age range.
For the model selected by lasso, the resulting model picks quite a small fraction of variables, with no interaction terms in the age_below_1 group, and one non-significant interaction terms in the age_above_1_below_5 group. Interaction term adv_mom:adv_rain is significant in the age_above_5 group while the other interaction term adv_mom:adv_sib is picked but not significant.
For the model selected by stepwise BIC, no interaction terms are selected in age_below_1 and age_above_1_below_5, which is the same result as in models selected by lasso. Meanwhile, adv_mom:adv_rain is picked and significant in age_above_5.
Interestingly, the cumulative model in the “age_below_1” group shows significance for the cumulative variable and the interaction term. There does not seem to be much worth noting in the other two groups.