In khellton/zibppca:

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

ZIBP-PCA: Zero-inflated Bivariate Poisson Principal Component analysis

This vignette describes how to use the R package zibppca to perform Zero-inflated Bivariate Poisson Principal Component analysis for Neuro-Psychiatric Inventory (NPI) items.

Installation

You can install the R package zibppca directly from GitHub, after installing the package devtools. The help description files can be accessed by using ?

#Install and load the devtools package
install.packages("devtools")
library(devtools)

# Install the zibppca package directly from GitHub
install_github("khellton/zibppca")

#Load the zibppca package
library(zibppca)

#Check the help file
?zibppca

Analyzing NPI items

First, load the NPI data either from SPSS or csv file and create a data.frame. Then the first 10 patients should look something like this:

library(foreign) # to import SPSS sav files
data_new <- read.spss('~/NR/PCANPI/Data/New_selection_nursing_home.sav', to.data.frame=TRUE)
data <- data.frame(data_new[(data_new$aar2004_2011==0),8:19])

get.names <- read.csv("~/NR/ModNPI/Data/NPILong.csv",colClasses = c('factor',rep('numeric',15)),dec =',',na.strings = '#NULL!', header = TRUE,sep = ';')
colnames(data) <- colnames(get.names)[4:15]
head(data, 10L)

The distribution of the different items can be visually inspected by histograms, or clustered box plots, for instance as seen below for Delusions and Aggression.

#Creating histograms
hist(data$Delusions_NPI1,main='Delusions',xlab = 'Domain score')
hist(data$Aggression_NPI3,main='Aggression',xlab = 'Domain score')

#Creating clustered box plots
par(mar=c(9,2,0,0))
boxplot(data,las=2,cex.axis=0.7)

Principal component analysis

The zibppca package has the following outline

library(zibppca)
#Run the main command 
result <- zibppca(data,scale. = TRUE)

### Loadings of the first 3 components
result$loadings[,1:3]

### Scree plot of eigenvalues
result$eigenvalues 
plot(result$eigenvalues,main = 'Scree plot',xlab='Eigenvalues')

### Plot patient scores 
plot(result$scores[,1],result$scores[,2],main = 'Score plot',xlab='1. component',ylab='2. component')