Markov Chain Monte Carlo (MCMC) to run Bayesian split population survival model with exchangeable frailties.
Returns a summary of a exchangeSPsurv object via summary.mcmc.
Print method for a exchangeSPsurv x.
Returns a plot of a exchangeSPsurv object via plot.mcmc.
exchangeSPsurv( duration, immune, Y0, LY, S, data, N, burn, thin, w = c(1, 1, 1), m = 10, ini.beta = 0, ini.gamma = 0, ini.W = 0, ini.V = 0, form = c("Weibull", "exponential", "loglog"), prop.varV, prop.varW, id_WV = unique(data[, S]) ) # S3 method for frailtySPsurv summary(object, parameter = character(), ...) # S3 method for frailtySPsurv print(x, ...) # S3 method for frailtySPsurv plot(x, parameter = character(), ...)
| duration | survival stage equation written in a formula of the form Y ~ X1 + X2 + ... where Y is duration until failure or censoring. |
|---|---|
| immune | split stage equation written in a formula of the form C ~ Z1 + Z2 + ... where C is a binary indicator of immunity. |
| Y0 | the elapsed time since inception until the beginning of time period (t-1). |
| LY | last observation year (coded as 1; 0 otherwise) due to censoring or failure. |
| S | spatial information (e.g. district ID) for each observation that matches the spatial matrix row/column information. |
| data | data.frame. |
| N | number of MCMC iterations. |
| burn | burn-in to be discarded. |
| thin | thinning to prevent from autocorrelation. |
| w | size of the slice in the slice sampling for (betas, gammas, rho). Write it as a vector. E.g. c(1,1,1). |
| m | limit on steps in the slice sampling. A vector of values for beta, gamma, rho. |
| ini.beta | initial value for the parameter vector beta. By default is 0. |
| ini.gamma | initial value for the parameter vector gamma. By default is 0. |
| ini.W | initial value for the parameter vector W. By default is 0. |
| ini.V | initial value for the parameter vector V. By default is 0. |
| form | type of parametric model (Weibull, Exponential, or Log-Logistic). |
| prop.varV | proposed variance for V in Metropolis-Hastings. |
| prop.varW | proposed variance for W in Metropolis-Hastings. |
| id_WV | vector of type character that modifies the colnames of W and V in the model’s result. By default is |
| object | an object of class |
| parameter | one of five parameters of the |
| ... | additional parameter. |
| x | an object of class |
exchangeSPsurv returns an object of class "frailtySPsurv".
An "exchangeSPsurv" object has the following elements:
matrix, numeric values of the posterior for each variable in the duration equation .
matrix, numeric values of the posterior for each variable in the immune equation.
vector, numeric values of rho.
vector, numeric values of lambda.
vector, numeric values of delta.
matrix, numeric values of the posterior for Ws.
matrix, numeric values of the posterior for Vs.
matrix of X's variables.
matrix of Z's variables.
vector of `Y'.
vector of `Y0'.
vector of `C'.
vector of `S'.
numeric initial values of beta.
numeric initial values of gamma.
numeric initial values of W.
numeric initial values of V.
character, type of distribution.
description for the model to be estimated.
# \donttest{ ## 1 walter <- spduration::add_duration(Walter_2015_JCR,"renewed_war", unitID = "ccode", tID = "year", freq = "year", ongoing = FALSE)#> Warning: Converting to 'Date' class with yyyy-06-30# add S walter <- spatial_SA(data = walter, var_ccode = "ccode", threshold = 800L) set.seed(123456) model <- exchangeSPsurv( duration = duration ~ fhcompor1 + lgdpl + comprehensive + victory + instabl + intensityln + ethfrac + unpko, immune = cured ~ fhcompor1 + lgdpl + victory, Y0 = 't.0', LY = 'lastyear', S = 'sp_id' , data = walter[[1]], N = 100, burn = 10, thin = 10, w = c(1,1,1), m = 10, form = "Weibull", prop.varV = 1e-05, prop.varW = 1e-05, ) print(model)#> Call: #> exchangeSPsurv(duration = duration ~ fhcompor1 + lgdpl + comprehensive + #> victory + instabl + intensityln + ethfrac + unpko, immune = cured ~ #> fhcompor1 + lgdpl + victory, Y0 = "t.0", LY = "lastyear", #> S = "sp_id", data = walter[[1]], N = 100, burn = 10, thin = 10, #> w = c(1, 1, 1), m = 10, form = "Weibull", prop.varV = 1e-05, #> prop.varW = 1e-05) #> #> #> Iterations = 1:9 #> Thinning interval = 1 #> Number of chains = 1 #> Sample size per chain = 9 #> #> Empirical mean and standard deviation for each variable, #> plus standard error of the mean: #> #> #> Duration equation: #> Mean SD Naive SE Time-series SE #> (Intercept) 1.37905266 0.66895063 0.22298354 0.22298354 #> fhcompor1 -1.05208201 0.41173391 0.13724464 0.13724464 #> lgdpl 0.09177184 0.09222484 0.03074161 0.03074161 #> comprehensive -0.62867094 0.59109153 0.19703051 0.18385024 #> victory 0.29303499 0.58880263 0.19626754 0.19626754 #> instabl 1.01191934 0.50979179 0.16993060 0.08099278 #> intensityln -0.02981777 0.08536595 0.02845532 0.02845532 #> ethfrac 0.21749487 0.71230113 0.23743371 0.09423386 #> unpko 0.42351151 0.40632993 0.13544331 0.13544331 #> #> Immune equation: #> Mean SD Naive SE Time-series SE #> (Intercept) 0.04137955 1.515039 0.5050131 0.5050131 #> fhcompor1 -0.62513545 2.151090 0.7170301 0.7170301 #> lgdpl -1.40500156 1.170416 0.3901386 0.3901386 #> victory 0.38827123 1.738965 0.5796550 0.5796550 #>#> #> Iterations = 1:9 #> Thinning interval = 1 #> Number of chains = 1 #> Sample size per chain = 9 #> #> 1. Empirical mean and standard deviation for each variable, #> plus standard error of the mean: #> #> Mean SD Naive SE Time-series SE #> (Intercept) 1.37905 0.66895 0.22298 0.22298 #> fhcompor1 -1.05208 0.41173 0.13724 0.13724 #> lgdpl 0.09177 0.09222 0.03074 0.03074 #> comprehensive -0.62867 0.59109 0.19703 0.18385 #> victory 0.29303 0.58880 0.19627 0.19627 #> instabl 1.01192 0.50979 0.16993 0.08099 #> intensityln -0.02982 0.08537 0.02846 0.02846 #> ethfrac 0.21749 0.71230 0.23743 0.09423 #> unpko 0.42351 0.40633 0.13544 0.13544 #> #> 2. Quantiles for each variable: #> #> 2.5% 25% 50% 75% 97.5% #> (Intercept) 0.40086 1.16042 1.272499 1.61321 2.45037 #> fhcompor1 -1.61784 -1.41128 -1.087209 -0.79473 -0.49926 #> lgdpl -0.02351 0.02864 0.079840 0.10775 0.24792 #> comprehensive -1.53418 -0.86159 -0.600034 -0.09388 0.06865 #> victory -0.56255 0.05853 0.218766 0.73028 1.05096 #> instabl 0.33299 0.74310 0.877161 1.38788 1.81246 #> intensityln -0.14131 -0.12867 -0.007338 0.01990 0.07484 #> ethfrac -0.56884 -0.45040 0.323703 0.64013 1.35275 #> unpko -0.25219 0.27497 0.380398 0.60618 0.98369 #># plot(model) ## 2 walter <- spduration::add_duration(Walter_2015_JCR,"renewed_war", unitID = "ccode", tID = "year", freq = "year", ongoing = FALSE)#> Warning: Converting to 'Date' class with yyyy-06-30walter$S <- rep(x = 1:length(unique(walter$ccode)), times = rle(walter$ccode)$lengths) country <- countrycode::countrycode(unique(walter$ccode),'gwn','iso3c') set.seed(123456) model <- exchangeSPsurv( duration = duration ~ fhcompor1 + lgdpl + comprehensive + victory + instabl + intensityln + ethfrac + unpko, immune = cured ~ fhcompor1 + lgdpl + victory, Y0 = 't.0', LY = 'lastyear', S = 'S' , data = walter, N = 100, burn = 10, thin = 10, w = c(1,1,1), m = 10, form = "loglog", prop.varV = 1e-05, prop.varW = 1e-05, id_WV = country ) print(model)#> Call: #> exchangeSPsurv(duration = duration ~ fhcompor1 + lgdpl + comprehensive + #> victory + instabl + intensityln + ethfrac + unpko, immune = cured ~ #> fhcompor1 + lgdpl + victory, Y0 = "t.0", LY = "lastyear", #> S = "S", data = walter, N = 100, burn = 10, thin = 10, w = c(1, #> 1, 1), m = 10, form = "loglog", prop.varV = 1e-05, prop.varW = 1e-05, #> id_WV = country) #> #> #> Iterations = 1:9 #> Thinning interval = 1 #> Number of chains = 1 #> Sample size per chain = 9 #> #> Empirical mean and standard deviation for each variable, #> plus standard error of the mean: #> #> #> Duration equation: #> Mean SD Naive SE Time-series SE #> (Intercept) 0.9768352 2.868855 0.9562849 2.4338006 #> fhcompor1 2.0920802 1.345390 0.4484634 0.7795867 #> lgdpl -3.6042466 2.072933 0.6909776 1.8029114 #> comprehensive 4.6996456 1.280050 0.4266832 1.0258938 #> victory 3.9664584 2.307090 0.7690301 1.4804038 #> instabl 1.9935549 2.001250 0.6670835 0.6670835 #> intensityln -2.1330003 1.285726 0.4285752 0.7461028 #> ethfrac -3.7590150 2.308131 0.7693771 1.5568944 #> unpko 4.0763098 1.558262 0.5194206 0.5194206 #> #> Immune equation: #> Mean SD Naive SE Time-series SE #> (Intercept) -0.7748578 1.1808708 0.3936236 0.3936236 #> fhcompor1 6.8669687 2.4097608 0.8032536 1.7212104 #> lgdpl -1.0312274 1.6514783 0.5504928 1.2353179 #> victory 1.2450188 0.9418921 0.3139640 0.3139640 #># }