A dataset containing items that measure Post-traumatic stress disorder symptoms (Armour et al. 2017). There are 20 variables (p) and 221 observations (n).
Estimation of Partial Correlation Matrix
The GGM can be estimated using estimate function with some important arguments:
type: could be continuous, binary, ordinal, or mixed
iter: Number of iterations (posterior samples; defaults to 5000).
prior_sd: Scale of the prior distribution, approximately the standard deviation of a beta distribution (defaults to 0.50).
Click to see the code
library(BGGM)library(tidyverse) # for ggplot2 and dplyrlibrary(psychonetrics) # the frequentist way for GGM# dataY <- ptsd[, 1:5] +1# add + 1 to make sure the ordinal variable start from '1'# ordinal modelfit <-estimate(Y, type ='ordinal', analytic =FALSE)fit_cont <-estimate(Y, type ='continuous', analytic =FALSE)#summary(fit)pcor_mat(fit)
freqEst <-pcor_to_df(getmatrix(mod_saturated_ord, 'omega')) freqEst_cont <-pcor_to_df(getmatrix(mod_saturated_cont, 'omega')) # transform partial corr. mat. of BGGM into df.BayesEst <-pcor_to_df(pcor_mat(fit))BayesEst_cont <-pcor_to_df(pcor_mat(fit_cont)) combEst <-rbind( freqEst |>mutate(Type ='psychonetrics_ord'), BayesEst |>mutate(Type ='BGGM_ord'), freqEst_cont |>mutate(Type ='psychonetrics_cont'), BayesEst_cont |>mutate(Type ='BGGM_cont')) |>mutate(Type =factor(Type, levels =c("psychonetrics_ord", "BGGM_ord", "psychonetrics_cont","BGGM_cont")))ggplot(combEst) +geom_col(aes(x =fct_rev(Label), y = Mean, fill = Type, col = Type), position =position_dodge()) +geom_text(aes(x =fct_rev(Label), y =ifelse(Mean>0, Mean+.03, Mean-.03), label =round(Mean, 3), group = Type), position =position_dodge(width = .9)) +labs(x ='') +coord_flip() +theme_classic()
We can find that:
Partial correlations with ordinal responses estimated by psychonetricdo not align well with partial correlations with ordinal responses estimated by BGGM
Partial correlations with continuous responses estimated by psychonetricalign well with partial correlations with continuous responses estimated by BGGM
Violating the continuous assumption of items responses usually underestimate the partial correlations (see BGGM_cont vs. BGGM_ord)
Regularization
It is easy to conduct regularization using Graph Selection using hypothesis testing
Click to see the code
fit_reg <- BGGM::select(fit, alternative ="greater")plot(fit_reg)
$plt
Hypothesis testing
Using BGGM, we can perform hypothesis testing like whether edge A is more strong than edge B, whether edge A is “significant” higher or smaller than 0, whether group A and group B have differences in partial correlations.
Testing Edge A vs. Edge B
Click to see the code
# example hypotheseshyp1 <-c("B2--B3 > B1--B4 > 0") # hypothesis 1hyp2 <-c("B1--B5 > B2--B5 > 0") # hypothesis 2hyp3 <-c("B2--B4 < 0") # hypothesis 3hyp4 <-c("B1--B3 < 0") # hypothesis 4## Posterior hypothesis probabilities.extract_php <-function(Y = Y, hypothesis = hyp, type ='ordinal') { x <-confirm(Y = Y, hypothesis = hypothesis, type = type) x_info <- x$info PHP_confirmatory <-as.matrix(x_info$PHP_confirmatory)rownames(PHP_confirmatory) <-c(x$hypothesis, "complement") PHP_confirmatory}extract_php(Y = Y, hypothesis = hyp1, type ='ordinal')
[,1]
B2--B3 > B1--B4 > 0 0.982
complement 0.018
Click to see the code
extract_php(Y = Y, hypothesis = hyp2, type ='ordinal')
[,1]
B1--B5 > B2--B5 > 0 0.855
complement 0.145
Click to see the code
extract_php(Y = Y, hypothesis = hyp3, type ='ordinal')
[,1]
B2--B4 < 0 0.686
complement 0.314
Click to see the code
extract_php(Y = Y, hypothesis = hyp4, type ='ordinal')
[,1]
B1--B3 < 0 0.272
complement 0.728
There are four hypothesis tested here:
The first expresses that the relation B2--B3 is larger than B1--B4 . In other words, that the partial correlation is larger for B2--B3. The p(H1|data) = 0.991 means the posterior hypothesis probability (PHP) is 0.991 which provides evidence for the hypothesis: the partial correlation for B2-B3 is larger;
The second expresses that the relation B1--B5 is larger than B2--B5 . The results of PHP suggests the hypothesis can be accepted;
The third expresses that the relation B2--B4 is larger than 0 . The results of PHP suggests the hypothesis can be accepted;
The fourth expresses that the relation B1--B3 is lower than 0 . The results of PHP suggests the hypothesis should be rejected given p(H_1) < p(complement);