R/stan_biglm.R
, R/stan_biglm.fit.R
stan_biglm.Rd
This is the same model as with stan_lm
but it utilizes the
output from biglm
in the biglm package in order to
proceed when the data is too large to fit in memory.
stan_biglm( biglm, xbar, ybar, s_y, ..., prior = R2(stop("'location' must be specified")), prior_intercept = NULL, prior_PD = FALSE, algorithm = c("sampling", "meanfield", "fullrank"), adapt_delta = NULL ) stan_biglm.fit( b, R, SSR, N, xbar, ybar, s_y, has_intercept = TRUE, ..., prior = R2(stop("'location' must be specified")), prior_intercept = NULL, prior_PD = FALSE, algorithm = c("sampling", "meanfield", "fullrank", "optimizing"), adapt_delta = NULL, importance_resampling = TRUE, keep_every = 1 )
biglm | The list output by |
---|---|
xbar | A numeric vector of column means in the implicit design matrix excluding the intercept for the observations included in the model. |
ybar | A numeric scalar indicating the mean of the outcome for the observations included in the model. |
s_y | A numeric scalar indicating the unbiased sample standard deviation of the outcome for the observations included in the model. |
... | Further arguments passed to the function in the rstan
package ( |
prior | Must be a call to |
prior_intercept | Either Note: If using a dense representation of the design matrix
---i.e., if the |
prior_PD | A logical scalar (defaulting to |
algorithm | A string (possibly abbreviated) indicating the
estimation approach to use. Can be |
adapt_delta | Only relevant if |
b | A numeric vector of OLS coefficients, excluding the intercept |
R | A square upper-triangular matrix from the QR decomposition of the design matrix, excluding the intercept |
SSR | A numeric scalar indicating the sum-of-squared residuals for OLS |
N | A integer scalar indicating the number of included observations |
has_intercept | A logical scalar indicating whether to add an intercept to the model when estimating it. |
importance_resampling | Logical scalar indicating whether to use
importance resampling when approximating the posterior distribution with
a multivariate normal around the posterior mode, which only applies
when |
keep_every | Positive integer, which defaults to 1, but can be higher
in order to thin the importance sampling realizations and also only
apples when |
The output of both stan_biglm
and stan_biglm.fit
is an
object of stanfit-class
rather than
stanreg-objects
, which is more limited and less convenient
but necessitated by the fact that stan_biglm
does not bring the full
design matrix into memory. Without the full design matrix,some of the
elements of a stanreg-objects
object cannot be calculated,
such as residuals. Thus, the functions in the rstanarm package that
input stanreg-objects
, such as
posterior_predict
cannot be used.
The stan_biglm
function is intended to be used in the same
circumstances as the biglm
function in the biglm
package but with an informative prior on the \(R^2\) of the regression.
Like biglm
, the memory required to estimate the model
depends largely on the number of predictors rather than the number of
observations. However, stan_biglm
and stan_biglm.fit
have
additional required arguments that are not necessary in
biglm
, namely xbar
, ybar
, and s_y
.
If any observations have any missing values on any of the predictors or the
outcome, such observations do not contribute to these statistics.
# create inputs ols <- lm(mpg ~ wt + qsec + am, data = mtcars, # all row are complete so ... na.action = na.exclude) # not necessary in this case b <- coef(ols)[-1] R <- qr.R(ols$qr)[-1,-1] SSR <- crossprod(ols$residuals)[1] not_NA <- !is.na(fitted(ols)) N <- sum(not_NA) xbar <- colMeans(mtcars[not_NA,c("wt", "qsec", "am")]) y <- mtcars$mpg[not_NA] ybar <- mean(y) s_y <- sd(y) post <- stan_biglm.fit(b, R, SSR, N, xbar, ybar, s_y, prior = R2(.75), # the next line is only to make the example go fast chains = 1, iter = 500, seed = 12345)#> #> SAMPLING FOR MODEL 'lm' NOW (CHAIN 1). #> Chain 1: #> Chain 1: Gradient evaluation took 1.9e-05 seconds #> Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.19 seconds. #> Chain 1: Adjust your expectations accordingly! #> Chain 1: #> Chain 1: #> Chain 1: Iteration: 1 / 500 [ 0%] (Warmup) #> Chain 1: Iteration: 50 / 500 [ 10%] (Warmup) #> Chain 1: Iteration: 100 / 500 [ 20%] (Warmup) #> Chain 1: Iteration: 150 / 500 [ 30%] (Warmup) #> Chain 1: Iteration: 200 / 500 [ 40%] (Warmup) #> Chain 1: Iteration: 250 / 500 [ 50%] (Warmup) #> Chain 1: Iteration: 251 / 500 [ 50%] (Sampling) #> Chain 1: Iteration: 300 / 500 [ 60%] (Sampling) #> Chain 1: Iteration: 350 / 500 [ 70%] (Sampling) #> Chain 1: Iteration: 400 / 500 [ 80%] (Sampling) #> Chain 1: Iteration: 450 / 500 [ 90%] (Sampling) #> Chain 1: Iteration: 500 / 500 [100%] (Sampling) #> Chain 1: #> Chain 1: Elapsed Time: 0.636559 seconds (Warm-up) #> Chain 1: 0.289041 seconds (Sampling) #> Chain 1: 0.9256 seconds (Total) #> Chain 1:#> Warning: There were 1 divergent transitions after warmup. Increasing adapt_delta above 0.99 may help. See #> http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup#> Warning: Examine the pairs() plot to diagnose sampling problems#> Warning: The largest R-hat is 1.06, indicating chains have not mixed. #> Running the chains for more iterations may help. See #> http://mc-stan.org/misc/warnings.html#r-hat#> 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 #> http://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 #> http://mc-stan.org/misc/warnings.html#tail-ess#> lm stan_lm #> wt -3.916504 -3.776792 #> qsec 1.225886 1.218030 #> am 2.935837 3.002232