spatFilt.bms             package:spatBMS             R Documentation

_B_a_y_e_s_i_a_n _M_o_d_e_l _S_a_m_p_l_i_n_g _a_n_d _A_v_e_r_a_g_i_n_g  _w_i_t_h _u_n_c_e_r_t_a_i_n _S_p_a_t_i_a_l _E_f_f_e_c_t_s

_D_e_s_c_r_i_p_t_i_o_n:

     Given data and prior information, this function samples all
     possible model combinations via MC3 and returns aggregate results.

_U_s_a_g_e:

     spatFilt.bms(X.data, WList,burn = 1000, iter = NA, nmodel = 100, mcmc = "bd",
       g = "UIP", mprior = "random", mprior.size = NA, user.int = TRUE, 
         start.value = NA, logfile = FALSE, logstep = 10000)

_A_r_g_u_m_e_n_t_s:

  X.data: a data frame or a matrix, with the dependent variable in the
          first column, followed by the covariates (alternatively,
          'X.data' can also be provided as a 'formula'). Note that
          'bms' automatically estimates a constant, therefore including
          constant terms is not necessary. 

   WList: a list object with length corresponding to the number of
          weight matrices you use ('length(WList)' must be '>1')). Each
          element of the list contains a matrix with the extracted
          eigenvectors, where the matrices do not have to have the same
          column dimension.  

    burn: The (positive integer) number of burn-in draws for the MC3
          sampler, defaults to 1000. (Not taken into account if
          mcmc="enumerate")

    iter: If mcmc is set to an MC3 sampler, then this is the number of
          iteration draws to be sampled (ex burn-ins), default 3000
          draws. 

  nmodel: the number of best models for which information is stored
          (default 500). Best models are used for convergence analysis
          between likelihoods and MCMC frequencies, as well as
          likelihood-based inference.
           Note that a very high value for 'nmodel' slows down the
          sampler significantly. Set nmodel=0 to speed up sampling (if
          best model information is not needed).

    mcmc: a character denoting the model sampler to be used.
           The MC3 sampler 'mcmc="bd"' corresponds to a birth/death
          MCMC algogrithm. 'mcmc="rev.jump"' enacts a reversible jump
          algorithm adding a "swap" step to the birth / death steps
          from "bd".
           Alternatively, the entire model space may be fully
          enumerated by setting 'mcmc="enumerate"' which will iterate
          all possible regressor combinations (Note: consider that this
          means 2^K iterations, where K is the number of covariates.)
           Cf. section 'Details' for more options. 

       g: the hyperparameter on Zellner's g-prior for the regression
          coefficients.
           'g="UIP"' corresponds to g=N, the number of observations
          (default);
           'g="BRIC"' corresponds to the benchmark prior suggested by
          Fernandez, Ley and Steel (2001), i.e  g=max(N, K^2), where K
          is the total number of covariates;
           'g="RIC"' sets g=K^2 and conforms to the risk inflation
          criterion by George and Foster (1994)
           'g="HQ"' sets g=log(N)^3 and asymptotically mimics the
          Hannan-Quinn criterion with C_{HQ}=3  (cf. Fernandez, Ley and
          Steel, 2001, p.395)


  mprior: a character denoting the model prior choice, defaulting to
          "random":
            'mprior="fixed"' denotes fixed common prior inclusion
          probabilities for each regressor as e.g. in Sala-i-Martin,
          Doppelhofer, and Miller(2004) -  for their fine-tuning, cf.
          'mprior.size'. Preferable to 'mcmc="random"' if strong prior
          information on model size exists;
           'mprior="random"' (default) triggers the 'random theta'
          prior by Ley and Steel (2008), who suggest a binomial-beta
          hyperprior on the a priori inclusion probability;
           'mprior="uniform"' employs the uniform model prior;
           'mprior="customk"' allows for custom model size priors (cf.
          'mprior.size');
           'mprior="pip"' allows for custom prior inclusion
          probabilities (cf. 'mprior.size');
           Note that the prior on models with more than N-3 regressors
          is automatically zero: these models will not be sampled.  

mprior.size: if 'mprior' is "fixed" or "random", 'mprior.size' is a
          scalar that denotes the prior expected value of the model
          size prior (default K/2).
           If 'mprior="customk"' then a custom model size prior can be
          provided as a K+1 vector detailing the priors from model size
          0 to K  (e.g. rep(1,K+1) for the uniform model prior);
           if 'mprior="pip"', then custom prior inclusion probabilities
          can be provided as a vector of size K, with elements in the
          interval (0,1) 

user.int: 'interactive mode': print out results to console after ending
          the routine and plots a chart (default TRUE). 

start.value: specifies the starting model of the iteration chain. For
          instance a specific model by the corresponding column indices
          (e.g. starting.model=numeric(K) starts from the null model
          including  solely a constant term) or 'start.value=c(3,6)'
          for a starting model only including covariates 3 and 6.
           If 'start.model' is set to an integer (e.g.
          'start.model=15') then that number of covariates (here: 15
          covariates) is randomly chosen and the starting model is
          identified by those regressors with an OLS t-statistic>0.2.
            The default value 'start.value=NA' corresponds to
          'start.value=min(ncol(X.data),nrow(X.data)-3)'. Note that
          'start.value=0' or 'start.value=NULL' starts from the null
          model.
           If 'mcmc="enumerate"' then 'start.value' is the index to
          start the iteration (default: 0, the null model) . Any number
          between 0 and K^2-1 is admissible. 

 logfile: setting 'logfile=TRUE' produces a logfile named '"test.log"'
          in your current working directory, in order to keep track of
          the sampling procedure. 'logfile' equal to some filepath
          (like 'logfile="subfolder/log.txt"') puts the logfile  into
          that specified position. (default: 'logfile=FALSE'). Note
          that 'logfile=""' implies log printouts on the console. 

 logstep: specifies at which number of posterior draws information is
          written to the log file; default: 10 000 iterations 

_D_e_t_a_i_l_s:

     Ad 'mcmc': 
      Interaction sampler: adding an ".int" to an MC3 sampler (e.g.
     "mcmc="bd.int") provides for special treatment of interaction
     terms. Interaction terms will only be sampled along with their
     component variables: In the colnumn names of X.data, interaction
     terms need to be  denominated by names consisting of the base
     terms separated by '#' (e.g. an interaction term of base variables
     '"A"', '"B"' and '"C"' needs column name '"A#B#C"'). Then variable
     '"A#B#C"' will only be included in a model if all of the component
     variables ("A", "B", and "C") are included. 

     The MC3 samplers "'bd'", "'rev.jump'", "'bd.int'" and
     "'rev.jump.int'", iterate away from a starting model by adding,
     droping or swapping (only in the case of rev.jump) covariates. 

     In an MCMC fashion, they thus randomly draw a candidate model and
     then move to it in case its marginal likelihood (marg.lik.) is
     superior to the marg.lik. of the current model. 

     In case the candidate's marg.lik is inferior, it is randomly
     accepted or rejected according to a probability formed by the
     ratio of candidate marg.lik over currrent marg.lik. Over time, the
     sampler should thus converge to a sensible distribution. For
     aggregate results based on these MC3 frequencies, the first few
     iterations are typically disregarded (the 'burn-ins').

_V_a_l_u_e:

     A list of classes 'bma' and 'spatFilt', that may be displayed
     using e.g. 'summary.bma' or 'coef.bma'. The list contains the
     following elements: 

    info: a list of aggregate statistics: 'iter' is the number of
          iterations, 'burn' the number of burn-ins.
           The following have to be divided by 'cumsumweights' to get
          posterior expected values: 'inccount' are the posterior
          inclusion probabilities, 'b1mo' and 'b2mo' the first and
          second moment of coefficients, 'add.otherstats' other
          statistics of interest (typically the moments of the
          shrinkage factor), 'msize' is the post. expected model size,
          'k.vec' the posterior model size distribution, 'pos.sign' the
          unconditional post. probability of positive coefficents,
          'corr.pmp' is the correlation between the best models' MCMC
          frequencies and their marg. likelihoods.
           'timed' is the time that was needed for MCMC sampling,
          'cons' is the posterior expected value of the constant. 'K'
          and 'N' are the maximum number of covariates and the sample
          size, respectively.

arguments: a list of the evaluated function arguments provided to 'bms'
          (see above)

  topmod: a 'topmod' object containing the best drawn models. see
          'topmod' for more details

start.pos: the positions of the starting model. If bmao is a'bma'
          object this corresponds to covariates
          bmao$reg.names[bmao$start.pos]. If bmao is a chain that
          resulted from several starting models (cf. 'c.bma', then
          'start.pos' is a list detailing all of them.

gprior.info: a list detailing information on the g-prior: 'gtype'
          corresponds to argument 'g' above, 'is.constant' is FALSE if
          'gtype' is either "hyper" or "EBL", 'return.g.stats'
          corresponds to argument 'g.stats' above, 'shrinkage.moments'
          contains the first and second moments of the shrinkage factor
          (only if 'return.g.stats==TRUE'), 'g' details the fixed g (if
          'is.constant==TRUE'), 'hyper.parameter' corresponds to the
          hyper-g parameter a as in Liang et al. (2008) 

mprior.info: a list detailing information on the model prior:
          'origargs' lists the original arguments to 'mprior' and
          'mprior.size' above; 'pmp(...)' is a function to calculate
          the prior model probability for a specific model; 'mp.mode'
          corresponds to argument 'mprior' above; 'mp.msize' denotes
          the prior mode size; 'mp.Kdist' is a (K+1) vector with the
          prior model size distribution from 0 to K

  X.data: data.frame or matrix: corresponds to argument 'X.data' above,
          possibly cleaned for NAs

reg.names: character vector: the covariate names to be used for X.data

spatFilt.call: the original call to the 'spatFilt.bms' function

wTopModels: a binary matrix containing the best nmodel models along an
          index (W-Index)  which W-matrix has been used.

   WList: corresponds to argument 'WList'.

  Wcount: a vector containing  the frequency with which the different
          weight matrices have been visited by the MCMC chain. Burnins
          are not counted.

_N_o_t_e:

     There are several ways to speed-up sampling: 'nmodel=10' saves
     only the ten best models, at most a marginal improvement.
     'nmodels=0' does not save the best (500) models,  however then
     posterior convergence and likelihood-based inference are not
     possible. 'g.stats=FALSE' saves some time by not retaining the
     shrinkage factors for the MC3 chain (and the best models).

_A_u_t_h_o_r(_s):

     Martin Feldkircher

_R_e_f_e_r_e_n_c_e_s:

     Crespo Cuaresma, Jesus and M. Feldkircher (2010): Spatial
     Filtering, Model Uncertainty and the Speed of Income Convergence
     in Europe. Working Papers 160, Oesterreichische Nationalbank.

     Feldkircher, M. and S. Zeugner (2009): Benchmark Priors Revisited:
     On Adaptive Shrinkage and the Supermodel Effect in Bayesian Model
     Averaging, IMF Working Paper 09/202.

     Fernandez, C. E. Ley and M. Steel (2001): Benchmark priors for
     Bayesian model averaging. Journal of Econometrics 100(2), 381-427 

     Ley, E. and M. Steel (2008): On the Effect of Prior Assumptions in
     Bayesian Model Averaging with Applications to Growth Regressions.
     working paper

     Tiefelsdorf, M. and D. Griffith (2007): Semiparametric filtering
     of spatial autocorrelation: the eigenvector approach. Environment
     and Planning A, 37, 1193-1221.

_S_e_e _A_l_s_o:

     'coef.bma', 'plotModelsize' and 'density.bma' for some operations
     on the resulting 'bma' object, 'c.bma' for integrating separate
     MC3 chains and splitting of sampling over several runs.

     Check <URL: http://bms.zeugner.eu> for additional help.

_E_x_a_m_p_l_e_s:

       library(BMS)
       data(dataBoston); data(WL.boston)
       #estimating a standard MC3 chain 
       bma1=spatFilt.bms(X.data=dataBoston,WList=WL.boston,burn=1e05,iter=1e05,
                   nmodel=100,mcmc="bd",g="bric",mprior="random",mprior.size=(ncol(dataBoston)-1)/2)
                   
       coef(bma1,exact=TRUE, std.coefs=TRUE) #standard coefficients based on exact likelihoods of the 100 best models
       
       #suppressing user-interactive output and  using a customized starting value (the null model)
       bma2 = spatFilt.bms(dataBoston,WList=WL.boston,burn=1000, iter=2000, nmodel=10, 
                           start.value=rep(0,ncol(dataBoston)),user.int=FALSE,)
       coef(bma2)
       
       #MC3 chain with the bric prior, saving only the 20 best models, and an alternative sampling procedure; putting a log entry to console every 1000th step
       bma3 = spatFilt.bms(dataBoston,burn=1000,WList=WL.boston, iter=5000, nmodel=20, g="bric", mcmc="rev.jump",logfile="",logstep=1000)
       image(bma3,cex=0.6) #showing the coefficient signs of the 20 best models and setting the font size of the labels in the plot to cex=0.6
       
       #using an interaction sampler for two interaction terms
       bma4 = spatFilt.bms(X.data=dataBoston,WList=WL.boston,burn=1000,iter=9000,start.value=0,mcmc="bd.int") 
       
       density(bma4,reg="NOX") # plot posterior density for covariate "NOX" with a grid
       grid()

