Nothing
### R code from vignette source 'examples.Rnw'
### Encoding: UTF-8
###################################################
### code chunk number 1: pima-data
###################################################
library(MASS)
pima <- rbind(Pima.tr, Pima.te)
pima$hasDiabetes <- as.numeric(pima$type == "Yes")
pima.nObs <- nrow(pima)
###################################################
### code chunk number 2: pima-setup
###################################################
library(hypergsplines)
modelData.pima <- with(pima,
glmModelData(y=hasDiabetes,
X=
cbind(npreg,
glu,
bp,
skin,
bmi,
ped,
age),
splineType="cubic",
nKnots=4L,
family=binomial))
###################################################
### code chunk number 3: pima-search-setup
###################################################
chainlength.pima <- 100
computation.pima <- getComputation(useOpenMP=FALSE,
higherOrderCorrection=FALSE)
###################################################
### code chunk number 4: pima-search
###################################################
set.seed(93)
time.pima <-
system.time(models.pima <-
stochSearch(modelData=modelData.pima,
modelPrior="dependent",
chainlength=chainlength.pima,
nModels=chainlength.pima,
computation=computation.pima))
###################################################
### code chunk number 5: pima-top-models
###################################################
head(models.pima$models)
map.pima <- models.pima$models[1, 1:7]
###################################################
### code chunk number 6: pima-incprobs
###################################################
round(models.pima$inclusionProbs,2)
###################################################
### code chunk number 7: pima-sampling
###################################################
mcmc.pima <- getMcmc(samples=500L,
burnin=100L,
step=1L,
nIwlsIterations=2L)
set.seed(634)
map.samples.pima <- glmGetSamples(config=map.pima,
modelData=modelData.pima,
mcmc=mcmc.pima,
computation=computation.pima)
###################################################
### code chunk number 8: pima-samples-structure
###################################################
str(map.samples.pima)
###################################################
### code chunk number 9: pima-plot-ex
###################################################
plotCurveEstimate(covName="age",
samples=map.samples.pima$samples,
modelData=modelData.pima)
###################################################
### code chunk number 10: pima-postprocess
###################################################
optim.map.pima <- postOptimize(modelData=modelData.pima,
modelConfig=map.pima,
computation=computation.pima)
optim.map.pima
###################################################
### code chunk number 11: pima-fit-samples
###################################################
fit.samples.pima <- getFitSamples(X=modelData.pima$origX[1:10,],
samples=map.samples.pima$samples,
modelData=modelData.pima)
obs.samples.pima <- plogis(fit.samples.pima)
###################################################
### code chunk number 12: pima-postpred-means
###################################################
rowMeans(obs.samples.pima)
###################################################
### code chunk number 13: pima-observations
###################################################
modelData.pima$Y[1:10]
###################################################
### code chunk number 14: pima-model-averaging
###################################################
average.samples.pima <-
with(models.pima,
getBmaSamples(config=models[1:10,],
logPostProbs=models$logMargLik[1:10] +
models$logPrior[1:10],
nSamples=500L,
modelData=modelData.pima,
mcmc=mcmc.pima,
computation=computation.pima))
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.