diff --git a/NEWS.md b/NEWS.md index 26f44162..68cdfe20 100644 --- a/NEWS.md +++ b/NEWS.md @@ -2,15 +2,15 @@ ## New Additions +### SAR models + The simultaneously-specified spatial autoregressive (SAR) model---referred to as the spatial error model (SEM) in the spatial econometrics literature---has been implemented. The SAR model can be applied directly to continuous data (as the likelihood function) or it can be used as prior model for spatially autocorrelated parameters. Details are provided on the documentation page for the `stan_sar` function. ## Minor changes -Previously, when getting fitted values from an auto-normal model (i.e., the CAR model with `family = auto_gaussian()`) the fitted values did not include the implicit spatial trend. Now, the `fitted.geostan_fit` method will return the fitted values with the implicit spatial trend; this is consistent with the behavior of `residuals.geostan_fit`, which has an option to `detrend` the residuals. This applies to the SAR and CAR auto-normal specifications. For details, see the documentation pages for `stan_car` and `stan_sar`. - -## Minor changes + - Previously, when getting fitted values from an auto-normal model (i.e., the CAR model with `family = auto_gaussian()`) the fitted values did not include the implicit spatial trend. Now, the `fitted.geostan_fit` method will return the fitted values with the implicit spatial trend; this is consistent with the behavior of `residuals.geostan_fit`, which has an option to `detrend` the residuals. This applies to the SAR and CAR auto-normal specifications. For details, see the documentation pages for `stan_car` and `stan_sar`. -The documentation for the models (`stan_glm`, `stan_car`, `stan_esf`, `stan_icar`, `stan_sar`) now uses Latex to typeset the model equations. + - The documentation for the models (`stan_glm`, `stan_car`, `stan_esf`, `stan_icar`, `stan_sar`) now uses Latex to typeset the model equations. # geostan 0.3.0 diff --git a/docs/articles/spatial-me-models.html b/docs/articles/spatial-me-models.html index 1c020751..ed95605d 100644 --- a/docs/articles/spatial-me-models.html +++ b/docs/articles/spatial-me-models.html @@ -221,9 +221,6 @@

Spatial ME model## Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable. ## Running the chains for more iterations may help. See ## https://mc-stan.org/misc/warnings.html#bulk-ess -
## Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
-## Running the chains for more iterations may help. See
-## https://mc-stan.org/misc/warnings.html#tail-ess

The CAR models in geostan can be highly efficient in terms of MCMC sampling, but the required number of iterations will depend on many characteristics of the model and data. Often the default iter = 2000 is more than sufficient (with cores = 4). To speed up sampling with multi-core processing, use cores = 4 (to sample 4 chains in parallel).

@@ -233,7 +230,7 @@

Visual diagnostics\(z_i\) an the modeled values \(x_i\): \[\delta_i = z_i - x_i.\] We want to look for any strong spatial pattern in these \(\delta_i\) values, because that would be an indication of a bias. However, the magnitude of the \(\delta_i\) value is important to consider—there may be a pattern, but if the amount of shrinkage is very small, that pattern may not matter. (The model returns \(M\) samples from the posterior distribution of each \(x_i\); or, indexing by MCMC sample \(x_i^m\) (\(m\) is an index, not an exponent). The reported \(\delta_i\) values is the MCMC mean \(\delta_i = \frac{1}{M} \sum_m z_i - x_i^m\).)

Two figures are provided to evaluate patterns in \(\delta_i\): first, a Moran scatter plot and, second, a map.

-
+
 me_diag(fit, 'insurance', georgia)

In this case, the results do not look too concerning insofar as there are no conspicuous patterns. However, a number of the credible intervals on the modeled values are large, which indicates that this data is not of high quality. The fact that some of the \(\delta_i\) values are substantial also points to low data quality.

@@ -243,54 +240,54 @@

Working with MCMC samples from

geostan consists of pre-compiled Stan models, and users can always access the Markov chain Monte Carlo (MCMC) samples returned by Stan. When extracted as a matrix of samples (as below), each row of the matrix represents a draw from the joint probability distribution for all model parameters, and each column consists of samples from the marginal distribution of each parameter.

The ME models return samples for \(x_i\) as well as the model parameters \(\mu\) (“mu_x_true”), \(\rho\) (“car_rho_x_true”), and \(\tau\) (“sigma_x_true”). We can access these using as.matrix (also as.array and as.data.frame).

-
+
 mu.x <- as.matrix(fit, pars = "mu_x_true")
 dim(mu.x)
## [1] 1300    1
-
+
 head(mu.x)
##           parameters
 ## iterations mu_x_true[1]
-##       [1,]     1.700061
-##       [2,]     1.866507
-##       [3,]     1.709768
-##       [4,]     1.721449
-##       [5,]     1.862860
-##       [6,]     1.771713
-
+##       [1,]     1.742830
+##       [2,]     1.745734
+##       [3,]     1.815951
+##       [4,]     1.830006
+##       [5,]     1.663939
+##       [6,]     1.677925
+
 mean(mu.x)
-
## [1] 1.784943
+
## [1] 1.786869

We can visualize these using plot or print a summary:

-
+
 print(fit$stanfit, pars = c("mu_x_true", "car_rho_x_true", "sigma_x_true"))
## Inference for Stan model: foundation.
 ## 4 chains, each with iter=650; warmup=325; thin=1; 
 ## post-warmup draws per chain=325, total post-warmup draws=1300.
 ## 
 ##                   mean se_mean   sd 2.5%  25%  50%  75% 97.5% n_eff Rhat
-## mu_x_true[1]      1.78    0.01 0.09 1.57 1.75 1.79 1.83  1.95   201 1.02
-## car_rho_x_true[1] 0.91    0.00 0.07 0.74 0.87 0.92 0.96  1.00   545 1.00
-## sigma_x_true[1]   0.48    0.00 0.04 0.42 0.46 0.48 0.51  0.56   496 1.00
+## mu_x_true[1]      1.79       0 0.07 1.63 1.75 1.79 1.83  1.93   517 1.01
+## car_rho_x_true[1] 0.91       0 0.06 0.75 0.87 0.92 0.95  0.99   658 1.00
+## sigma_x_true[1]   0.48       0 0.04 0.42 0.46 0.48 0.51  0.57   631 1.00
 ## 
-## Samples were drawn using NUTS(diag_e) at Mon Sep 19 15:02:35 2022.
+## Samples were drawn using NUTS(diag_e) at Mon Sep 19 18:58:29 2022.
 ## 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).

To extract samples from the joint probability distribution for \(\boldsymbol x\), use the generic parameter name “x_true”:

-
+
 x <- as.matrix(fit, pars = "x_true")
 dim(x)
## [1] 1300  159

If we wanted to calculate the mean of each of these marginal distributions (one for every \(x_i\)), we could use apply with MARGIN = 2 to summarize by column:

-
+
 x.mu <- apply(x, 2, mean)
 head(x.mu)
## x_insurance[1] x_insurance[2] x_insurance[3] x_insurance[4] x_insurance[5] 
-##      0.8364901      0.7999472      0.8518636      0.8523471      0.9179334 
+##      0.8365794      0.7994134      0.8519202      0.8521496      0.9178745 
 ## x_insurance[6] 
-##      0.8702361
+## 0.8701874

The vector x.mu contains estimates (posterior means) for \(x_i\). We might want to use these to plot the residuals or fitted values against the predictor:

-
+
 rs <- resid(fit)$mean
 plot(x.mu, rs)
 abline(h = 0)
@@ -301,7 +298,7 @@

Working with MCMC samples from

Non-spatial ME models

If the ME list doesn’t have a slot with car_parts, geostan will automatically use a non-spatial Student’s t model instead of the CAR model: \[ p(\boldsymbol x | \mathcal{M}) = Student(\boldsymbol x | \nu, \mu \boldsymbol 1, \sigma) \]

-
+
 ME_nsp <- prep_me_data(
   se = data.frame(insurance = georgia$insurance.se),
   logit = TRUE
diff --git a/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-14-1.png b/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-14-1.png
index fbc07313..0ed525a4 100644
Binary files a/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-14-1.png and b/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-14-1.png differ
diff --git a/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-9-1.png b/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-9-1.png
index 93344481..1283abd8 100644
Binary files a/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-9-1.png and b/docs/articles/spatial-me-models_files/figure-html/unnamed-chunk-9-1.png differ
diff --git a/docs/news/index.html b/docs/news/index.html
index 4bf11db6..62b5a8e0 100644
--- a/docs/news/index.html
+++ b/docs/news/index.html
@@ -60,16 +60,16 @@ 

Changelog

New Additions

+
+

SAR models

The simultaneously-specified spatial autoregressive (SAR) model—referred to as the spatial error model (SEM) in the spatial econometrics literature—has been implemented. The SAR model can be applied directly to continuous data (as the likelihood function) or it can be used as prior model for spatially autocorrelated parameters. Details are provided on the documentation page for the stan_sar function.

-
-

Minor changes

-

Previously, when getting fitted values from an auto-normal model (i.e., the CAR model with family = auto_gaussian()) the fitted values did not include the implicit spatial trend. Now, the fitted.geostan_fit method will return the fitted values with the implicit spatial trend; this is consistent with the behavior of residuals.geostan_fit, which has an option to detrend the residuals. This applies to the SAR and CAR auto-normal specifications. For details, see the documentation pages for stan_car and stan_sar.

-

Minor changes

-

The documentation for the models (stan_glm, stan_car, stan_esf, stan_icar, stan_sar) now uses Latex to typeset the model equations.

-
+

Minor changes

+
  • Previously, when getting fitted values from an auto-normal model (i.e., the CAR model with family = auto_gaussian()) the fitted values did not include the implicit spatial trend. Now, the fitted.geostan_fit method will return the fitted values with the implicit spatial trend; this is consistent with the behavior of residuals.geostan_fit, which has an option to detrend the residuals. This applies to the SAR and CAR auto-normal specifications. For details, see the documentation pages for stan_car and stan_sar.

  • +
  • The documentation for the models (stan_glm, stan_car, stan_esf, stan_icar, stan_sar) now uses Latex to typeset the model equations.

  • +
diff --git a/docs/pkgdown.yml b/docs/pkgdown.yml index 177aa5cc..58fe3924 100644 --- a/docs/pkgdown.yml +++ b/docs/pkgdown.yml @@ -4,7 +4,7 @@ pkgdown_sha: ~ articles: measuring-sa: measuring-sa.html spatial-me-models: spatial-me-models.html -last_built: 2022-09-19T20:48Z +last_built: 2022-09-20T00:06Z urls: reference: https://connordonegan.github.io/geostan/reference article: https://connordonegan.github.io/geostan/articles diff --git a/docs/reference/geostan_fit.html b/docs/reference/geostan_fit.html index 16d02c47..a48ada5f 100644 --- a/docs/reference/geostan_fit.html +++ b/docs/reference/geostan_fit.html @@ -72,7 +72,7 @@

geostan_fit methods

) # S3 method for geostan_fit -plot(x, pars, plotfun = "hist", fill = "steelblue4", ...) +plot(x, pars, plotfun = "hist", fill = "steelblue4", ...) # S3 method for geostan_fit as.matrix(x, ...) @@ -184,7 +184,7 @@

Examples

# print and plot results print(fit) -plot(fit) +plot(fit) # residuals r = resid(fit) @@ -202,7 +202,7 @@

Examples

# posterior predictive distribution yrep <- posterior_predict(fit, S = 65) -plot(density(yrep[1,]), col = "gray30") +plot(density(yrep[1,]), col = "gray30") for (i in 2:nrow(yrep)) lines(density(yrep[i,]), col = "gray30") lines(density(georgia$deaths.male), col = "darkred", lwd = 2) @@ -217,7 +217,7 @@

Examples

) p <- predict(fit, newdata, type = "response") -plot(newdata$income, p$mean * 1e3, +plot(newdata$income, p$mean * 1e3, type = 'l', main = "Deaths per 1,000", ylab = NA,