---
title: "Getting Started with bgms"
output: rmarkdown::html_vignette
vignette: >
  %\VignetteIndexEntry{Getting Started with bgms}
  %\VignetteEngine{knitr::rmarkdown}
  %\VignetteEncoding{UTF-8}
---

```{r, include = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>",
  fig.width = 6,
  fig.height = 4
)
```

# Introduction

The **bgms** package implements Bayesian methods for analyzing graphical
models of binary and ordinal variables. It estimates main effects
(category thresholds) and pairwise interactions in an ordinal Markov
random field (MRF), with optional Bayesian edge selection via
spike–and–slab priors. The package provides two main entry points:

- `bgm()` for one-sample designs (single network),
- `bgmCompare()` for independent-sample designs (group comparisons).

This vignette walks through the basic workflow: fitting a model,
summarizing posterior output, and visualizing results.


# Wenchuan dataset

The dataset `Wenchuan` contains responses from survivors of the 2008 Wenchuan
earthquake on posttraumatic stress items. Here, we analyze a subset of the
first five items as a demonstration.


```{r}
library(bgms)

# Analyse a subset of the Wenchuan dataset
?Wenchuan
data = Wenchuan[, 1:5]
head(data)
```

# Fitting a model

The main entry point is `bgm()` for single-group models and 
`bgmCompare()` for multiple-group comparisons.


```{r, include=FALSE}
fit <- readRDS(system.file("extdata", "fit_5items.rds", package = "bgms"))
```

```{r, eval=FALSE}
fit = bgm(data, seed = 1234)
```
Note: During fitting, progress bars are shown in interactive sessions.
In this vignette, they are suppressed for clarity. Sampling can take a while; 
the progress bars usually help track progress.

# Posterior summaries

```{r}
summary(fit)
```

You can also access posterior means or inclusion probabilities directly:

```{r}
coef(fit)
```

# Network plot

To visualize the network structure, we threshold the posterior inclusion probabilities at 0.5 and plot the resulting adjacency matrix.

```{r, fig.width= 7, fig.height= 7}
library(qgraph)

median_probability_network = coef(fit)$pairwise
median_probability_network[coef(fit)$indicator < 0.5] = 0.0

qgraph(median_probability_network, 
       theme = "TeamFortress", 
       maximum = 1,
       fade = FALSE,
       color = c("#f0ae0e"), vsize = 10, repulsion = .9, 
       label.cex = 1, label.scale = "FALSE", 
       labels = colnames(data))
```

# Next steps

- For comparing groups, see `?bgmCompare` or the *Model Comparison* vignette.
- For diagnostics and convergence checks, see the *Diagnostics* vignette.