# Simulate two historical datasets
n <- 100
P <- 4
time1 <- round(rexp(n, rate=0.5),1)
event1 <- rep(1,n)
X1 <- matrix(rbinom(n*P,prob=0.5,size=1), ncol=P)
S1 <- c(rep(1,n/2),rep(2,n/2))
time2 <- round(rexp(n, rate=0.7),1)
event2 <- rep(1,n)
X2 <- matrix(rbinom(n*P,prob=0.5,size=1), ncol=P)
S2 <- c(rep(1,n/2),rep(2,n/2))
historical <- list(list(time=time1, event=event1, X=X1, S=S1),
list(time=time2, event=event2, X=X2, S=S2))
# We choose three intervals for the first stratum and two intervals for the second stratum
n.intervals <- c(3,2)
change.points <- list(c(1,2), 2)
# Get samples from the approximate normalized power prior for beta
nMC <- 100 # nMC should be larger in practice
nBI <- 50
prior.beta <- approximate.prior.beta(historical, n.intervals, change.points=change.points,
prior.a0.shape1=c(1,1), prior.a0.shape2=c(1,1),
nMC=nMC, nBI=nBI)
prior_beta_mu=colMeans(prior.beta)
prior_beta_sigma=cov(prior.beta)
# Aprroximate the discrete sames with a single multivariate normal with weight one.
# The user can use any mixture of multivariate normal distributions as an
# approximation for the normalized power prior for beta.
prior.beta.mvn <- list(list(prior_beta_mu, prior_beta_sigma, 1))
# prior.beta.mvn is a parameter for phm.random.a0() and power.phm.random.a0()
Run the code above in your browser using DataLab