Description Usage Arguments Details Value Author(s) Examples
Interfaces to spBayes
functions that can be used
in a pipeline implemented by magrittr
.
1 2 3 |
data |
data frame, tibble, list, ... |
... |
Other arguments passed to the corresponding interfaced function. |
Interfaces call their corresponding interfaced function.
Object returned by interfaced function.
Roberto Bertolusso
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 | ## Not run:
library(intubate)
library(magrittr)
library(spBayes)
## NOTE: Just interfacing 3 functions as a proof of concept. However,
## there are so many things declared and used, that I am not sure if
## everything could be passed on a pipeline (perhaps with intuBags down the line).
## I may get back to it later.
## ntbt_bayesGeostatExact: Simple Bayesian spatial linear model with fixed semivariogram parameters
data(FORMGMT.dat)
n <- nrow(FORMGMT.dat)
p <- 5 ##an intercept an four covariates
n.samples <- 50
phi <- 0.0012
coords <- cbind(FORMGMT.dat$Longi, FORMGMT.dat$Lat)
coords <- coords*(pi/180)*6378
beta.prior.mean <- rep(0, times=p)
beta.prior.precision <- matrix(0, nrow=p, ncol=p)
alpha <- 1/1.5
sigma.sq.prior.shape <- 2.0
sigma.sq.prior.rate <- 10.0
## Original function to interface
bayesGeostatExact(Y ~ X1 + X2 + X3 + X4, data = FORMGMT.dat,
n.samples = n.samples,
beta.prior.mean = beta.prior.mean,
beta.prior.precision = beta.prior.precision,
coords = coords, phi = phi, alpha = alpha,
sigma.sq.prior.shape = sigma.sq.prior.shape,
sigma.sq.prior.rate = sigma.sq.prior.rate)
## The interface puts data as first parameter
ntbt_bayesGeostatExact(FORMGMT.dat, Y ~ X1 + X2 + X3 + X4,
n.samples = n.samples,
beta.prior.mean = beta.prior.mean,
beta.prior.precision = beta.prior.precision,
coords = coords, phi = phi, alpha = alpha,
sigma.sq.prior.shape = sigma.sq.prior.shape,
sigma.sq.prior.rate = sigma.sq.prior.rate)
## so it can be used easily in a pipeline.
FORMGMT.dat %>%
ntbt_bayesGeostatExact(Y ~ X1 + X2 + X3 + X4,
n.samples = n.samples,
beta.prior.mean = beta.prior.mean,
beta.prior.precision = beta.prior.precision,
coords = coords, phi = phi, alpha = alpha,
sigma.sq.prior.shape = sigma.sq.prior.shape,
sigma.sq.prior.rate = sigma.sq.prior.rate)
## ntbt_bayesLMConjugate: Simple Bayesian linear model via the Normal/inverse-Gamma conjugate
n <- nrow(FORMGMT.dat)
p <- 7 ##an intercept and six covariates
n.samples <- 500
## Below we demonstrate the conjugate function in the special case
## with improper priors. The results are the same as for the above,
## up to MC error.
beta.prior.mean <- rep(0, times=p)
beta.prior.precision <- matrix(0, nrow=p, ncol=p)
prior.shape <- -p/2
prior.rate <- 0
## Original function to interface
bayesLMConjugate(Y ~ X1 + X2 + X3 + X4 + X5 + X6, data = FORMGMT.dat,
n.samples, beta.prior.mean,
beta.prior.precision,
prior.shape, prior.rate)
## The interface puts data as first parameter
ntbt_bayesLMConjugate(FORMGMT.dat, Y ~ X1 + X2 + X3 + X4 + X5 + X6,
n.samples, beta.prior.mean,
beta.prior.precision,
prior.shape, prior.rate)
## so it can be used easily in a pipeline.
FORMGMT.dat %>%
ntbt_bayesLMConjugate(Y ~ X1 + X2 + X3 + X4 + X5 + X6,
n.samples, beta.prior.mean,
beta.prior.precision,
prior.shape, prior.rate)
## ntbt_spDynLM: Function for fitting univariate Bayesian dynamic
## space-time regression models
data("NETemp.dat")
ne.temp <- NETemp.dat
set.seed(1)
##take a chunk of New England
ne.temp <- ne.temp[ne.temp[,"UTMX"] > 5500000 & ne.temp[,"UTMY"] > 3000000,]
##subset first 2 years (Jan 2000 - Dec. 2002)
y.t <- ne.temp[,4:27]
N.t <- ncol(y.t) ##number of months
n <- nrow(y.t) ##number of observation per months
##add some missing observations to illistrate prediction
miss <- sample(1:N.t, 10)
holdout.station.id <- 5
y.t.holdout <- y.t[holdout.station.id, miss]
y.t[holdout.station.id, miss] <- NA
##scale to km
coords <- as.matrix(ne.temp[,c("UTMX", "UTMY")]/1000)
max.d <- max(iDist(coords))
##set starting and priors
p <- 2 #number of regression parameters in each month
starting <- list("beta"=rep(0,N.t*p), "phi"=rep(3/(0.5*max.d), N.t),
"sigma.sq"=rep(2,N.t), "tau.sq"=rep(1, N.t),
"sigma.eta"=diag(rep(0.01, p)))
tuning <- list("phi"=rep(5, N.t))
priors <- list("beta.0.Norm"=list(rep(0,p), diag(1000,p)),
"phi.Unif"=list(rep(3/(0.9*max.d), N.t), rep(3/(0.05*max.d), N.t)),
"sigma.sq.IG"=list(rep(2,N.t), rep(10,N.t)),
"tau.sq.IG"=list(rep(2,N.t), rep(5,N.t)),
"sigma.eta.IW"=list(2, diag(0.001,p)))
##make symbolic model formula statement for each month
mods <- lapply(paste(colnames(y.t),'elev',sep='~'), as.formula)
n.samples <- 10 # in original example it is 2000.
## Original function to interface
spDynLM(mods, data=cbind(y.t,ne.temp[,"elev",drop=FALSE]), coords=coords,
starting=starting, tuning=tuning, priors=priors, get.fitted =TRUE,
cov.model="exponential", n.samples=n.samples, n.report=25)
## The interface puts data as first parameter
ntbt_spDynLM(cbind(y.t,ne.temp[,"elev",drop=FALSE]), mods, coords=coords,
starting=starting, tuning=tuning, priors=priors, get.fitted =TRUE,
cov.model="exponential", n.samples=n.samples, n.report=25)
## so it can be used easily in a pipeline.
cbind(y.t,ne.temp[,"elev",drop=FALSE]) %>%
ntbt_spDynLM(mods, coords=coords,
starting=starting, tuning=tuning, priors=priors, get.fitted =TRUE,
cov.model="exponential", n.samples=n.samples, n.report=25)
## End(Not run)
|
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.