Technical Supplement
- Probability Integral Transformation
Let \(Y\) be a continuous random variable with df \(F\) and quantile function \(F^{-1}\) and define \(F(Y)=Z\) then \(P(Z \le z) = P(F(Y) \le z)= P(F^{-1}(F(Y)) \le F^{-1}(z)) = P(Y \le F^{-1}(z))=F(F^{-1}(z)) =z\) which is recognized as the df of a standard uniform random variable. So \(Z \sim Unif(0,1)\).
- Formula for C-volume with arbitrary \(d\)
where the sum is over the set of \(2^d\) vectors of size \(d\) with elements 0 or 1. For \(d=2\), the set is \(\{(i_1=0,i_2=0),(i_1=0,i_2=1),(i_1=1,i_2=0),(i_1=1,i_2=1)\}\) and the C-volume is \(\Delta_{(\mathbf{a},\mathbf{b}]}C=C(b_1,b_2)-C(b_1,a_2)-C(a_1,b_2)+C(a_1,a_2)\). For \(d=3\) the set is \(\{(i_1=0,i_2=0,i_3=0),(i_1=0,i_2=0,i_3=1),(i_1=0,i_2=1,i_3=0),(i_1=0,i_2=1,i_3=1),(i_1=1,i_2=0,i_3=0),(i_1=1,i_2=0,i_3=1),(i_1=1,i_2=1,i_3=0),(i_1=1,i_2=1,i_3=1)\}\) and the C-volume is \(\Delta_{(\mathbf{a},\mathbf{b}]}C= C(b_1,b_2,b_3)- C(b_1,b_2,a_3)- C(b_1,a_2,b_3)+ C(b_1,a_2,a_3)- C(a_1,b_2,b_3)+ C(a_1,b_2,a_3)+ C(a_1,a_2,b_3)- C(a_1,a_2,a_3)\).
- Likelihood for discrete margins
- Likelihood for mixed margins
Where \(c_{\theta,S_1}\) is the copula density for the continuous variables (with the discrete margins marginalized out) and \(C_{\theta,S_2|S_1}\) is the copula for the discrete margins conditional on the continuous margins. The first line of equation (9) is similar to the log-likelihood equation for continuous variables in (5), but uses a copula density where the discrete margins have been averaged out. The second line of equation (9) is similar to the finite differencing for discrete variables in (8) with the modification that the copula is conditioned on the continuous margins.
In addition to the discrete and mixed cases shown above, Joe [13] also provides log-likelihoods for margins with right-censoring and missing-at-random values. Gunawan et al. [70] develop extensions for the case when an individual margin has both a discrete and continuous component.
- Mixed continuous/binary margins likelihood from Benefit-Risk example
- Benefit-Risk MCMC model diagnostics
## marginal model diagnostics
#http://mc-stan.org/bayesplot/
print(fit0a, probs=c(0.025, 0.975))
## Inference for Stan model: 968550e09f4560174cecbd127bc5bb9b.
## 4 chains, each with iter=3000; warmup=1000; thin=1;
## post-warmup draws per chain=2000, total post-warmup draws=8000.
##
## mean se_mean sd 2.5% 97.5% n_eff Rhat
## beta1[1] -153.00 0.11 9.77 -171.97 -133.92 8000 1
## beta1[2] -49.34 0.11 10.05 -68.81 -29.60 8000 1
## sigma 98.41 0.06 5.01 89.22 108.95 7838 1
## lp__ -1017.61 0.02 1.25 -1020.87 -1016.20 4033 1
##
## Samples were drawn using NUTS(diag_e) at Mon Feb 25 17:33:03 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).
color_scheme_set("viridis")
stan_trace(fit0a)
print(fit0b, probs=c(0.025, 0.975))
## Inference for Stan model: 14df10564cc86d04063140f36084b07d.
## 4 chains, each with iter=3000; warmup=1000; thin=1;
## post-warmup draws per chain=2000, total post-warmup draws=8000.
##
## mean se_mean sd 2.5% 97.5% n_eff Rhat
## beta2[1] -0.96 0.00 0.15 -1.26 -0.67 6695 1
## beta2[2] -0.28 0.00 0.13 -0.53 -0.03 7374 1
## p[1] 0.17 0.00 0.04 0.10 0.25 6547 1
## p[2] 0.39 0.00 0.05 0.30 0.49 7374 1
## lp__ -113.47 0.02 1.02 -116.14 -112.49 3453 1
##
## Samples were drawn using NUTS(diag_e) at Mon Feb 25 17:33:58 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).
stan_trace(fit0b)
## joint model diagnostics
print(br_fit, probs=c(0.025, 0.975))
## Inference for Stan model: 35135559ebf8c381d1a6ad94ea441ba4.
## 4 chains, each with iter=3000; warmup=1000; thin=1;
## post-warmup draws per chain=2000, total post-warmup draws=8000.
##
## mean se_mean sd 2.5% 97.5% n_eff Rhat
## beta1[1] -152.97 0.11 10.17 -172.84 -133.05 8000 1
## beta1[2] -49.60 0.13 9.88 -68.88 -30.33 5976 1
## s[1] 99.68 0.08 7.28 86.58 115.08 8000 1
## s[2] 97.69 0.08 6.88 85.07 112.14 8000 1
## beta2[1] -0.95 0.00 0.15 -1.25 -0.67 8000 1
## beta2[2] -0.30 0.00 0.13 -0.55 -0.05 5464 1
## omega[1] 0.02 0.00 0.14 -0.26 0.29 8000 1
## omega[2] 0.77 0.00 0.06 0.63 0.88 7168 1
## mu[1] -152.97 0.11 10.17 -172.84 -133.05 8000 1
## mu[2] -49.60 0.13 9.88 -68.88 -30.33 5976 1
## p[1] 0.17 0.00 0.04 0.11 0.25 8000 1
## p[2] 0.38 0.00 0.05 0.29 0.48 5486 1
## theta[1] 0.02 0.00 0.14 -0.26 0.29 8000 1
## theta[2] 0.77 0.00 0.06 0.63 0.88 7168 1
## rho[1] 0.01 0.00 0.09 -0.17 0.19 8000 1
## rho[2] 0.61 0.00 0.05 0.49 0.69 7066 1
## lp__ -1294.09 0.03 2.02 -1298.88 -1291.16 3835 1
##
## Samples were drawn using NUTS(diag_e) at Mon Feb 25 17:35:10 2019.
## For each parameter, n_eff is a crude measure of effective sample size,
## and Rhat is the potential scale reduction factor on split chains (at
## convergence, Rhat=1).
br_posterior <- extract(br_fit, inc_warmup = TRUE, permuted = FALSE)
color_scheme_set("mix-blue-red")
mcmc_trace(br_posterior, pars = c("mu[1]", "mu[2]", "p[1]", "p[2]", "rho[1]", "rho[2]",
"theta[1]", "theta[2]"), n_warmup = 1000,
facet_args = list(nrow = 4, labeller = label_parsed))
br_posterior2 <- extract(br_fit, inc_warmup = FALSE, permuted = FALSE)
mcmc_areas(br_posterior2, pars = c("mu[1]", "mu[2]"),
prob = 0.95, prob_outer = 0.99, point_est = "mean")
mcmc_areas(br_posterior2, pars = c("p[1]", "p[2]", "rho[1]", "rho[2]",
"theta[1]", "theta[2]"),
prob = 0.95, prob_outer = 0.99, point_est = "mean")
References
13. Joe H. Dependence modeling with copulas. Boca Raton: CRC Press, Taylor & Francis Group; 2015.
70. Gunawan D, Khaled MA, Kohn R. Mixed marginal copula modeling. Journal of Business & Economic Statistics [Internet] 2018 [cited 2018 Oct 12];1–11. Available from: https://www.tandfonline.com/doi/full/10.1080/07350015.2018.1469998