Effect plot in R

 library(effects)

x1 <- predictorEffect("x6", model)

plot(x1)

x2 <- predictorEffect("x2", model) 

plot(x2, main="model COVID 19",xlab="X2", ylab="x")  #single plot

plot(allEffects(model)) #all plot

model summary in R

 # check that model is good fit or not

with(model, cbind(res.deviance = deviance, df = df.residual,

                   p = pchisq(deviance, df.residual, lower.tail=FALSE)))

## odds ratios and 95% CI ***********

exp(cbind(OR = coef(model), confint(model)))

# MORE SUMMARIES ###########################################

anova(model)            # Coefficients w/inferential tests

coef(model)             # Coefficients 

hist(coef(model))

confint(model)          # CI for coefficients

hist(confint(model))

resid(model)            # Residuals case-by-case

hist(residuals(model),main="model COVID 19" ) # Histogram of residuals

plot(residuals(model), main="model COVID 19" )

logLik(model)

BIC(model)

PseudoR2(model)

predict(model)

hist(predict(model))

#peseudo r square 

model$null.deviance

model$deviance

modelChi <- model$null.deviance - model$deviance

pseudo <- modelChi / model$null.deviance

pseudo

# Compute the pseudo p-value

Chidf <- model$df.null - model$df.residual

modelChi <- model$null.deviance - model$deviance

1 - pchisq(modelChi, Chidf)

#RSS(residual sum of square)

RSS <- c(crossprod(model$residuals))

RSS

#Mean square error

MSE <- RSS / length(model$residuals)

MSE

#Root MSE

RMSE <- sqrt(MSE)

RMSE

#Pearson estimated residual variance

sig2 <- RSS / model$df.residual

sig2

Poisson lognormal regression model in R

 library(PLNmodels)

library(ggplot2)

library(corrplot)

Y=read.csv(choose.files())

X=read.csv(choose.files())

c=list(Y,X)

names(c) <- c("output", "input")

d <- prepare_data(c$output, c$input)

model <- PLN(Abundance ~ 

             LANE_WIDTH               +

             LENGTH                    +

             LOG_AVG_AADT                 +

             LOG_PAVEMENT_CONDITION         +

             HIGH_FREQ_TRANSIT           +

             eightyfivepercent_SPEED_ADJ   +

             P_SPEEDING_ADJ               +

             SL                         +

             NL, data = d)

print(model)      #variational lower bound of the ICL

coef(model)       #mu -vectors of means of the latent variable

sigma(model)      #covariance matrix of the latent variable

vcov(model)       #Variance-Covariance Matrix

fitted(model)

standard_error(model)

plot(fitted(model))

barplot(fitted(model))

#2 tailed z test

z<- coef(model)/standard_error(model)

p<- (1-pnorm(abs(z),0,1))*2

p        #pvalue

p value for any statistical model using R

 #2 tailed z test

z<- coef(model)/standard_error(model)

p<- (1-pnorm(abs(z),0,1))*2

p        #pvalue