Nothing
## -----------------------------------------------------------------------------
cov <- rnorm(500, 0.2)
cov2 <- rnorm(500, 4)
sp.2 <- ceiling(rnorm(500, 1) + (cov * 2))
sp.2[sp.2 < 0] <- 0
poiss.dat <- data.frame(sp.1 = ceiling(rnorm(500, 4) + (cov2 * 15.5) + (sp.2 * -0.15)),
sp.2 = sp.2, sp.3 = ceiling((sp.2 * 1.5) + rnorm(500, 0.1)))
poiss.dat[poiss.dat < 0] <- 0
poiss.dat$cov <- cov
poiss.dat$cov2 <- cov2
## -----------------------------------------------------------------------------
apply(poiss.dat[, c(1:3)], 2, range)
## -----------------------------------------------------------------------------
sp.1.vs.sp.2 <- coef(glm(sp.1 ~ sp.2, family = 'poisson', data = poiss.dat))[2]
sp.2.vs.sp.1 <- coef(glm(sp.2 ~ sp.1, family = 'poisson', data = poiss.dat))[2]
mean.coef <- mean(sp.1.vs.sp.2, sp.2.vs.sp.1)
# Exponentiate, due to the log link
mean.coef <- exp(mean.coef)
# Make plots of predicted vs observed
# First, predict sp.1 (the common species) abundances
plot(poiss.dat$sp.1, poiss.dat$sp.2 * mean.coef,
xlab = 'Observed',
ylab = 'Predicted')
# Now predict sp.2 (the more rare species) abundances
plot(poiss.dat$sp.2, poiss.dat$sp.1 * mean.coef,
xlab = 'Observed',
ylab = 'Predicted')
## -----------------------------------------------------------------------------
library(MRFcov)
poiss.crf <- MRFcov(data = poiss.dat, n_nodes = 3, family = 'poisson')
poiss.crf$poiss_sc_factors
## -----------------------------------------------------------------------------
sd(poiss.dat[,1]) * poiss.crf$direct_coefs$cov2[1]
## -----------------------------------------------------------------------------
poiss.preds <- predict_MRF(data = poiss.dat, MRF_mod = poiss.crf,
n_cores = 1)
plot(poiss.dat$sp.1, poiss.preds[,1],
xlab = 'Observed',
ylab = 'Predicted')
plot(poiss.dat$sp.2, poiss.preds[,2],
xlab = 'Observed',
ylab = 'Predicted')
## -----------------------------------------------------------------------------
poiss.cv <- cv_MRF_diag(data = poiss.dat, n_nodes = 3,
n_folds = 5,
n_cores = 1, family = 'poisson',
compare_null = TRUE, plot = FALSE)
# CRF (with covariates) model deviance
range(poiss.cv$Deviance[poiss.cv$model == 'CRF'])
# MRF (no covariates) model deviance
range(poiss.cv$Deviance[poiss.cv$model != 'CRF'])
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.