In this tutorial, we are going to demonstrate the relcs package. Using relcs, you are able to easily specify random effects latent change score models that allow proportional change (or, self-feedback) to vary between persons.
First, download the package from GitHub and attach it to your workspace. Here we also load a couple of other libraries that will make our life less miserable.
require("relcs") require("ggplot2") library("ggthemes") library(tidyr) set.seed(543231)
Let's start with simulating some data from a latent change score model. Here, we assume proportional change (mean=.2) and individual differences in the proportional change (variance=.1). We assume no additional linear slope over time. Residual error variance is set to 1. There is some individual differences at the intercept (mu=5, variance=1).
N <- 100 num.obs <- 5 simulated.data <- relcs::simulateDataFromRELCS(N,num.obs,selffeedback.mean = .2, selffeedback.variance = .05,has.slope = FALSE, residualerrorvariance = 1, interceptmu = 5, interceptvariance = 1) summary(simulated.data)
Let's first create a long version of our simulated data set, which we will use for plotting shortly. The tidyr package is a great helper for this task.
simulated.data <- cbind(ID=rep(1:N), simulated.data) long.simulated.data <- gather(data = simulated.data,X,value,-ID) long.simulated.data$ID <- as.factor(long.simulated.data$ID)
We'll use the amazing ggplot2 library to plot the data
subdata <- long.simulated.data ggplot(data=subdata,aes(x=X,y=value,color=ID,group=ID))+geom_line()+ guides(colour=FALSE)+ theme_fivethirtyeight()
Now, we will fit a RELCS model. Yes, this is really all you need to do:
fitted.model <- fitRELCS(data=simulated.data)
fitted.model <- fitRELCS(data=simulated.data, iter=200, warmup=40)
Yes, it is as simple as that! Let's have a look at the parameter point estimates (here, the means of the posterior distribution):
getEstimates(fitted.model)
Now, let's plot our estimates using relcs's in-built plotting functions.
est <- getEstimates(fitted.model) plot(est)
Now, realistically, your data will have some missing values here and there. Let's simulate a matrix that reflects the pattern of missing and non-missing values and remove some of our observations according to it. We set the probability of missing data to 10% here.
prob.miss <- 0.1 miss.mat <- ifelse(matrix(data = rbinom(n = prod(dim(simulated.data)),size = 1,prob=prob.miss),nrow=nrow(simulated.data)),NA,1) missing.data <- miss.mat * simulated.data
And again, we simply pass this data object to fitRELCS and it will handle the missing data appropriately:
require(relcs) fitted.model <- fitRELCS(data=missing.data, type="stan")
require(relcs) fitted.model <- fitRELCS(data=missing.data, type="stan")
Let's again check the estimates using print and plot:
est <- getEstimates(fitted.model) print(est) plot(est)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.