Last data update: 2014.03.03

R: Extract a Model from a bma Object
as.zlmR Documentation

Extract a Model from a bma Object

Description

Extracts a model out of a bma object's saved models and converts it to a zlm linear model

Usage

as.zlm(bmao, model = 1)

Arguments

bmao

A bma object, e.g. resulting from a call to bms

model

The model index, in one of the following forms:
An integer, denoting the rank of the model (1 for best, 2 for second-best, ...)
A numeric or logical vector of length K describing which covariates are contained in the model
A hexcode character describing which covariates are contained in the model

Details

A bma object stores several 'best' models it encounters (cf. argument nmodel in bms). as.zlm extracts a single model and converts it to an object of class zlm, which represents a linear model estimated under Zellner's g prior.
The utility model.frame allows to transfrom a zlm model into an OLS model of class lm.

Value

a list of class zlm

Author(s)

Stefan Zeugner

See Also

bms for creating bma objects, zlm for creating zlm objects, topmodels.bma and pmp.bma for displaying the topmodels in a bma object

Check http://bms.zeugner.eu for additional help.

Examples

data(datafls)

mm=bms(datafls[,1:6],mcmc="enumeration") # do a small BMA chain
topmodels.bma(mm)[,1:5] #display the best 5 models

m2a=as.zlm(mm,4) #extract the fourth best model
summary(m2a)

# Bayesian Model Selection:
# transform the best model into an OLS model:
lm(model.frame(as.zlm(mm)))

# extract the model only containing the 5th regressor
m2b=as.zlm(mm,c(0,0,0,0,1)) 

# extract the model only containing the 5th regressor in hexcode
print(bin2hex(c(0,0,0,0,1)))
m2c=as.zlm(mm,"01")




Results


R version 3.3.1 (2016-06-21) -- "Bug in Your Hair"
Copyright (C) 2016 The R Foundation for Statistical Computing
Platform: x86_64-pc-linux-gnu (64-bit)

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> library(BMS)
> png(filename="/home/ddbj/snapshot/RGM3/R_CC/result/BMS/as.zlm.Rd_%03d_medium.png", width=480, height=480)
> ### Name: as.zlm
> ### Title: Extract a Model from a bma Object
> ### Aliases: as.zlm
> ### Keywords: models
> 
> ### ** Examples
> 
> data(datafls)
> 
> mm=bms(datafls[,1:6],mcmc="enumeration") # do a small BMA chain
               PIP     Post Mean      Post SD Cond.Pos.Sign Idx
Abslat   0.5696540  0.0001690225 0.0001796401     1.0000000   1
WarDummy 0.4531863 -0.0040773389 0.0053271202     0.0000000   5
Spanish  0.4117183 -0.0045103513 0.0064646326     0.0000000   2
French   0.3816369 -0.0048899949 0.0074391511     0.0000000   3
Brit     0.1684193  0.0002939201 0.0026550967     0.5542162   4

Mean no. regressors               Draws             Burnins                Time 
           "1.9846"                "32"                 "0"  "0.005727291 secs" 
 No. models visited      Modelspace 2^K           % visited         % Topmodels 
               "32"                "32"               "100"               "100" 
           Corr PMP            No. Obs.         Model Prior             g-Prior 
               "NA"                "72"      "random / 2.5"               "UIP" 
    Shrinkage-Stats 
        "Av=0.9863" 

Time difference of 0.005727291 secs
> topmodels.bma(mm)[,1:5] #display the best 5 models
                  10        00         0d         01         1d
Abslat      1.000000 0.0000000 0.00000000 0.00000000 1.00000000
Spanish     0.000000 0.0000000 1.00000000 0.00000000 1.00000000
French      0.000000 0.0000000 1.00000000 0.00000000 1.00000000
Brit        0.000000 0.0000000 0.00000000 0.00000000 0.00000000
WarDummy    0.000000 0.0000000 1.00000000 1.00000000 1.00000000
PMP (Exact) 0.186113 0.1199236 0.07438732 0.05812197 0.05726452
PMP (MCMC)  0.186113 0.1199236 0.07438732 0.05812197 0.05726452
> 
> m2a=as.zlm(mm,4) #extract the fourth best model
> summary(m2a)
Coefficients
               Exp.Val.     St.Dev.
(Intercept)  0.02495729          NA
WarDummy    -0.01049913 0.004234111

 Log Marginal Likelihood:
133.7997
 g-Prior: UIP 
Shrinkage Factor: 0.986
> 
> # Bayesian Model Selection:
> # transform the best model into an OLS model:
> lm(model.frame(as.zlm(mm)))

Call:
lm(formula = model.frame(as.zlm(mm)))

Coefficients:
(Intercept)       Abslat  
  0.0115755    0.0003557  

> 
> # extract the model only containing the 5th regressor
> m2b=as.zlm(mm,c(0,0,0,0,1)) 
> 
> # extract the model only containing the 5th regressor in hexcode
> print(bin2hex(c(0,0,0,0,1)))
[1] "01"
> m2c=as.zlm(mm,"01")
> 
> 
> 
> 
> 
> 
> 
> 
> 
> dev.off()
null device 
          1 
>