PPD: Function to predict the joint probability involving two event...

Description Usage Arguments Details Value Author(s) References See Also Examples

View source: R/PPD.R

Description

PPD is a function to predict the joint probability involving two event times in Bayesian illness-death models.

Usage

1
PPD(fit, x1, x2, x3, t1, t2)

Arguments

fit

an object of class Bayes_HReg. Currently, the function is available for PEM illness-death models.

x1

a vector of covariates for h_1 with which to predict.

x2

a vector of covariates for h_2 with which to predict.

x3

a vector of covariates for h_3 with which to predict.

t1

time to non-terminal event for which the joint probability is calculated.

t2

time to terminal event for which the joint probability is calculated.

Details

Using the posterior predictive density, given (x_1, x_2, x_3), one can predict any joint probability involving the two event times such as P(T_1<t_1, T_2<t_2| x_1, x_2, x_3) for 0<t_1≤ t_2 and P(T_1=∞, T_2<t_2| x_1, x_2, x_3) for t_2>0.

Value

F_u

Predicted P(T_1≤ t_1, T_2≤ t_2| x_1, x_2, x_3) in the upper wedge of the support of (T_1, T_2).

F_l

Predicted P(T_1=∞, T_2≤ t_2| x_1, x_2, x_3) in the lower wedge of the support of (t1, t2).

Author(s)

Kyu Ha Lee and Sebastien Haneuse
Maintainer: Kyu Ha Lee <klee15239@gmail.com>

References

Lee, K. H., Haneuse, S., Schrag, D., and Dominici, F. (2015), Bayesian semiparametric analysis of semicompeting risks data: investigating hospital readmission after a pancreatic cancer diagnosis, Journal of the Royal Statistical Society: Series C, 64, 2, 253-273.

See Also

BayesID_HReg

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
## Not run:    
# loading a data set
data(scrData)
id=scrData$cluster

form <- Formula(time1 + event1 | time2 + event2 ~ x1 + x2 | x1 + x2 | x1 + x2)

#####################
## Hyperparameters ##
#####################

## Subject-specific frailty variance component
##  - prior parameters for 1/theta
##
theta.ab <- c(0.7, 0.7)

## PEM baseline hazard function
##
PEM.ab1 <- c(0.7, 0.7) # prior parameters for 1/sigma_1^2
PEM.ab2 <- c(0.7, 0.7) # prior parameters for 1/sigma_2^2
PEM.ab3 <- c(0.7, 0.7) # prior parameters for 1/sigma_3^2
##
PEM.alpha1 <- 10 # prior parameters for K1
PEM.alpha2 <- 10 # prior parameters for K2
PEM.alpha3 <- 10 # prior parameters for K3

##
hyperParams <- list(theta=theta.ab,
                   PEM=list(PEM.ab1=PEM.ab1, PEM.ab2=PEM.ab2, PEM.ab3=PEM.ab3,
                    PEM.alpha1=PEM.alpha1, PEM.alpha2=PEM.alpha2,
                    PEM.alpha3=PEM.alpha3))
                    
###################
## MCMC SETTINGS ##
###################

## Setting for the overall run
##
numReps    <- 2000
thin       <- 10
burninPerc <- 0.5

## Settings for storage
##
nGam_save <- 0

## Tuning parameters for specific updates
##
##  - those common to all models
mhProp_theta_var  <- 0.05
##
## - those specific to the Weibull specification of the baseline hazard functions
mhProp_alphag_var <- c(0.01, 0.01, 0.01)
##
## - those specific to the PEM specification of the baseline hazard functions
Cg        <- c(0.2, 0.2, 0.2)
delPertg  <- c(0.5, 0.5, 0.5)
rj.scheme <- 1
Kg_max    <- c(50, 50, 50)
sg_max    <- c(max(scrData$time1[scrData$event1 == 1]),
               max(scrData$time2[scrData$event1 == 0 & scrData$event2 == 1]),
               max(scrData$time2[scrData$event1 == 1 & scrData$event2 == 1]))

time_lambda1 <- seq(1, sg_max[1], 1)
time_lambda2 <- seq(1, sg_max[2], 1)
time_lambda3 <- seq(1, sg_max[3], 1)               

##
mcmc.PEM <- list(run=list(numReps=numReps, thin=thin, burninPerc=burninPerc),
                storage=list(nGam_save=nGam_save),
                tuning=list(mhProp_theta_var=mhProp_theta_var,
                Cg=Cg, delPertg=delPertg,
                rj.scheme=rj.scheme, Kg_max=Kg_max,
                time_lambda1=time_lambda1, time_lambda2=time_lambda2,
                time_lambda3=time_lambda3))
    
##						
myModel <- c("semi-Markov", "PEM")
myPath  <- "Output/02-Results-PEM/"

startValues      <- initiate.startValues_HReg(form, scrData, model=myModel, nChain=2)

##
fit_PEM <- BayesID_HReg(form, scrData, id=NULL, model=myModel,
                 hyperParams, startValues, mcmc.PEM, path=myPath)
				
PPD(fit_PEM, x1=c(1,1), x2=c(1,1), x3=c(1,1), t1=3, t2=6)
		

## End(Not run)

SemiCompRisks documentation built on Feb. 3, 2021, 5:06 p.m.