The R package simPreg simulates pregnancy data using
hazard rates derived from Norwegian registry data as default inputs.
Pregnancies can end in spontaneous or non-spontaneous live birth from
week 22, or in late miscarriage/stillbirth from week 12. Because the
distribution of early miscarriages (before week 12) is uncertain, these
are omitted from the simulation. A time-varying exposure is also
simulated, with default hazard rates derived from the distribution of
COVID-19 vaccination timing during pregnancy in Norway in 2020–2023.
The function simPregProp() calculates pregnancy outcome
proportions across gestational age and exposure timing within a
time-to-event framework.
library(simPreg)
df_prop <- simPregProp()
head(df_prop)
#> GA ExpGA Outcome Prop
#> 1 154 NA spontaneous_live_birth 1.391126e-05
#> 2 155 NA spontaneous_live_birth 4.633503e-06
#> 3 156 NA spontaneous_live_birth 6.943707e-06
#> 4 157 NA spontaneous_live_birth 6.937910e-06
#> 5 158 NA spontaneous_live_birth 1.848539e-05
#> 6 159 NA spontaneous_live_birth 4.617685e-06The function simPregSamp() generates a sample of
user-specified size n from the proportions produced by
simPregProp(). The logical argument expand
determines whether the sampled frequencies are returned as grouped
counts or expanded to one row per pregnancy.
Example using expand = FALSE:
set.seed(33)
df_samp_grp <- simPregSamp(
df = df_prop,
n = 100000,
expand = FALSE
)
head(df_samp_grp)
#> GA ExpGA Outcome Prop Freq
#> 1 154 NA spontaneous_live_birth 1.391126e-05 1
#> 2 157 NA spontaneous_live_birth 6.937910e-06 2
#> 3 158 NA spontaneous_live_birth 1.848539e-05 3
#> 4 160 NA spontaneous_live_birth 9.227821e-06 1
#> 5 161 NA spontaneous_live_birth 5.532308e-05 4
#> 6 164 NA spontaneous_live_birth 3.219793e-05 4Example using expand = TRUE:
set.seed(33)
df_samp_one <- simPregSamp(
df = df_prop,
n = 100000,
expand = TRUE
)
head(df_samp_one)
#> GA ExpGA Outcome Prop Freq
#> 1 154 NA spontaneous_live_birth 1.391126e-05 1
#> 2 157 NA spontaneous_live_birth 6.937910e-06 1
#> 3 157 NA spontaneous_live_birth 6.937910e-06 1
#> 4 158 NA spontaneous_live_birth 1.848539e-05 1
#> 5 158 NA spontaneous_live_birth 1.848539e-05 1
#> 6 158 NA spontaneous_live_birth 1.848539e-05 1
table(df_samp_one$Outcome)
#>
#> late_miscarriage_stillbirth nonspontaneous_live_birth
#> 484 29122
#> spontaneous_live_birth
#> 70394All hazard rate vectors can be replaced by user-specified vectors, provided that they have the same length. The default vectors have length 301 and are derived by survival analysis of real-world data.
By default, no effect of exposure on any pregnancy outcome is used in the simulations. Users can specify hazard ratios that are constant throughout gestation (proportional hazards) or vary over gestation.
Specifying an increased (proportional) hazard of late miscarriage/stillbirth following exposure:
# HR = 5
df_prop_sb <- simPregProp(hr.late.miscarriage.stillbirth = rep(5, 301))
set.seed(34)
df_samp_sb <- simPregSamp(
df = df_prop_sb,
n = 100000,
expand = TRUE
)
# Calculate proportions among exposed and unexposed pregnancies
prop_sb_exposed <- mean(
df_samp_sb$Outcome[!is.na(df_samp_sb$ExpGA)] ==
"late_miscarriage_stillbirth"
)
prop_sb_unexposed <- mean(
df_samp_sb$Outcome[is.na(df_samp_sb$ExpGA)] ==
"late_miscarriage_stillbirth"
)
c(
exposed = round(prop_sb_exposed, 3),
unexposed = round(prop_sb_unexposed, 3)
)
#> exposed unexposed
#> 0.016 0.005A non-null hazard ratio for preterm birth following exposure can be specified by setting the HR to a value different from 1 during the preterm period.
# HR = 2 from gestational day 154 through day 258, and 1 otherwise
hr_preterm <- c(rep(1, 153), rep(2, 105), rep(1, 43))
df_prop_pt <- simPregProp(
hr.spont.livebirth = hr_preterm,
hr.nonspont.livebirth = hr_preterm
)
set.seed(35)
df_samp_pt <- simPregSamp(
df = df_prop_pt,
n = 100000,
expand = TRUE
)
# Calculate proportions among exposed and unexposed pregnancies
prop_pt_exposed <- mean(
df_samp_pt$GA[!is.na(df_samp_pt$ExpGA)] >= 154 &
df_samp_pt$GA[!is.na(df_samp_pt$ExpGA)] <= 258
)
prop_pt_unexposed <- mean(
df_samp_pt$GA[is.na(df_samp_pt$ExpGA)] >= 154 &
df_samp_pt$GA[is.na(df_samp_pt$ExpGA)] <= 258
)
c(
exposed = round(prop_pt_exposed, 3),
unexposed = round(prop_pt_unexposed, 3)
)
#> exposed unexposed
#> 0.084 0.049The simulated data can be analyzed with a Cox proportional hazards
model from the survival package, treating exposure as
time-varying. Pregnancies are considered unexposed before exposure and
exposed thereafter, with follow-up represented in start–stop format.
# Exclude late miscarriage/stillbirth
df_samp_pt <- subset(
df_samp_pt,
Outcome != "late_miscarriage_stillbirth"
)
# Define preterm birth
df_samp_pt$preterm <- as.integer(df_samp_pt$GA <= 258)
# Start follow-up at day 153 to allow events from day 154 through day 258,
# and censor at day 258 otherwise
df_samp_pt$tstart <- 153
df_samp_pt$tstop <- pmin(df_samp_pt$GA, 258)
# Assign unique identifiers
df_samp_pt$id <- seq_len(nrow(df_samp_pt))
# Create start-stop data with exposure as a time-varying covariate
df_samp_pt_tv <- survival::tmerge(
data1 = df_samp_pt,
data2 = df_samp_pt,
id = id,
tstart = tstart,
tstop = tstop,
preterm = event(tstop, preterm),
exposed = tdc(ExpGA)
)
# Fit Cox proportional hazards model
cox_pt <- survival::coxph(
survival::Surv(tstart, tstop, preterm) ~ exposed,
data = df_samp_pt_tv
)
# Display estimated HR and 95% CI
c(
HR = round(exp(coef(cox_pt)), 3),
lower_95 = round(exp(confint(cox_pt)[1]), 3),
upper_95 = round(exp(confint(cox_pt)[2]), 3)
)
#> HR.exposed lower_95 upper_95
#> 1.989 1.874 2.112This example demonstrates how simPreg can be used to
simulate pregnancy data with time-varying exposures and specified
exposure effects. The simulated data can be analyzed using standard
time-to-event methods.