Loading Libraries

Loading libraries used in this script.

library(psych)     #describing the data
library(plyr)      #data manipulation
library(tidyr)     #tidy data
library(dplyr)     #data manipulation
library(ggplot2)   #data visualization
library(entropy)   #entropy calculation
library(nFactors)  #factor analysis

Loading Data

Our example makes use of the person-level and interaction-level (EMA-type) AMIB data files. We make use of person-level personality variables and a small set of interaction-level variables.

Loading person-level file and subsetting to variables of interest

#set filepath for data file
filepath <- "https://raw.githubusercontent.com/The-Change-Lab/collaborations/main/AMIB/AMIB_persons.csv"
#read in the .csv file using the url() function
AMIB_persons <- read.csv(file=url(filepath), header=TRUE)

#subsetting to variables of interest
AMIB_persons <- AMIB_persons %>%
  select(id, bfi_o, bfi_c, bfi_e, bfi_a, bfi_n)

Loading day-level file (T = 21) and subsetting to variables of interest.

#set filepath for data file
filepath <- "https://raw.githubusercontent.com/The-Change-Lab/collaborations/refs/heads/main/AMIB/AMIB_dailyP2.csv"
#read in the .csv file using the url() function
AMIB_daily <- read.csv(file=url(filepath),header=TRUE)

Note, while the AMIB person level data has 629 participants, the day level data has 30.

Calculating univariate IIV metrics

istats_univariate <- AMIB_daily %>%
  dplyr::group_by(id) %>%   #specify the unit we want to summarise by, in this case per person
  dplyr::summarize(icount_posaff = sum(!is.na(posaff)),  #count of observations
            imean_posaff  = mean(posaff, na.rm=TRUE), #imean (continuous)
            isd_posaff    = sd(posaff, na.rm=TRUE),   #isd   (continuous)
            iskew_posaff  = skew(posaff, na.rm=TRUE), #iskew (continuous)
            ikurt_posaff  = kurtosi(posaff, na.rm=TRUE), # ikurtosis (continuous)
            ientropy_stress= entropy(table(stress, useNA="no")))  #(categorical)

Look at univariate IIV metrics data

round(head(istats_univariate),2)

## # A tibble: 6 × 7
##      id icount_posaff imean_posaff isd_posaff iskew_posaff ikurt_posaff
##   <dbl>         <dbl>        <dbl>      <dbl>        <dbl>        <dbl>
## 1   203            22         4.42       0.84        -0.7          0.13
## 2   204            22         3.1        1.17        -0.58        -1.19
## 3   205            22         4.32       0.68        -0.73        -0.53
## 4   208            22         3.92       0.67         0.11        -0.86
## 5   211            22         2.81       0.85         0.7         -0.15
## 6   214            15         4.16       1.01        -0.46        -0.82
## # ℹ 1 more variable: ientropy_stress <dbl>

Describe the univariate IIV metrics

describe(istats_univariate[-1])

##                 vars  n  mean   sd median trimmed  mad   min   max range  skew
## icount_posaff      1 30 20.90 2.58  22.00   21.58 0.00 14.00 22.00  8.00 -2.02
## imean_posaff       2 30  3.78 0.71   3.91    3.76 0.56  2.61  5.62  3.01  0.29
## isd_posaff         3 30  0.88 0.26   0.87    0.87 0.26  0.32  1.52  1.20  0.17
## iskew_posaff       4 30 -0.15 0.46  -0.24   -0.19 0.46 -0.74  0.97  1.71  0.69
## ikurt_posaff       5 30 -0.79 0.50  -0.88   -0.84 0.47 -1.43  0.58  2.01  0.94
## ientropy_stress    6 30  1.47 0.28   1.50    1.48 0.28  0.84  1.90  1.06 -0.39
##                 kurtosis   se
## icount_posaff       2.28 0.47
## imean_posaff       -0.17 0.13
## isd_posaff         -0.19 0.05
## iskew_posaff       -0.45 0.08
## ikurt_posaff        0.01 0.09
## ientropy_stress    -0.77 0.05

Examine correlations

round(cor(istats_univariate[ ,-1], use="pairwise.complete.obs"),2)

##                 icount_posaff imean_posaff isd_posaff iskew_posaff ikurt_posaff
## icount_posaff            1.00         0.08      -0.32         0.16         0.24
## imean_posaff             0.08         1.00      -0.45        -0.38         0.14
## isd_posaff              -0.32        -0.45       1.00         0.21        -0.07
## iskew_posaff             0.16        -0.38       0.21         1.00         0.29
## ikurt_posaff             0.24         0.14      -0.07         0.29         1.00
## ientropy_stress         -0.06        -0.24       0.29        -0.21        -0.11
##                 ientropy_stress
## icount_posaff             -0.06
## imean_posaff              -0.24
## isd_posaff                 0.29
## iskew_posaff              -0.21
## ikurt_posaff              -0.11
## ientropy_stress            1.00

Merge the univariate IIV metrics with other person-level data

AMIB_univariate <- istats_univariate %>%
  left_join(AMIB_persons, by="id")

Examine relation between positive affect variability and neuroticism

Look at association between isd_posaff and bfi_n:

AMIB_univariate %>%
  ggplot(aes(x=bfi_n, y=isd_posaff)) +
  geom_point(colour="gray40") +
  geom_smooth(aes(group=1), method=lm, se=TRUE, fullrange=TRUE, lty=1,
              linewidth=2, color="#9A23DE") +
  scale_y_continuous(breaks=c(0.4, 0.8, 1.2, 1.6), limits=c(0.32,1.6)) +
  xlab("Neuroticism") + ylab("Positive Affect Variability") +
  theme_classic() +
  theme(axis.title=element_text(size=16),
        axis.text=element_text(size=12),
        plot.title=element_text(size=16, hjust=.5, face="bold")) +
  ggtitle("Between-Person Association Univariate IIV")

And the formal model test (cross-sectional regression).

reg_univariate <- lm(isd_posaff ~ bfi_n, 
                     data=AMIB_univariate,
                     na.action = na.exclude) 
summary(reg_univariate)

## 
## Call:
## lm(formula = isd_posaff ~ bfi_n, data = AMIB_univariate, na.action = na.exclude)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.49597 -0.14662 -0.00169  0.18189  0.49217 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.59383    0.15329   3.874 0.000589 ***
## bfi_n        0.08654    0.04499   1.923 0.064654 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.252 on 28 degrees of freedom
## Multiple R-squared:  0.1167, Adjusted R-squared:  0.08516 
## F-statistic: 3.699 on 1 and 28 DF,  p-value: 0.06465

• For the prototypical person, positive affect variability is expected to be 0.59
• Individuals who are higher on neuroticism tend to have more variability in their positive affect. For every unit increase in neuroticism, characteristic positive affect variability is expected to increase by 0.09 (however, this association is not statistically significantly different from zero p = 0.06).

Calculate bivariate IVV metrics

Calculate intraindiviudal correlation (iCor) and intraindividual covariance (iCov) and intraindividual regression (iReg):

istats_bivariate <- AMIB_daily %>%
  dplyr::group_by(id) %>%
  dplyr::summarise(icor_affect = cor(x=posaff, y=negaff, 
                              use="pairwise.complete.obs", method="pearson"),
            icov_affect = cov(x=posaff, y=negaff, 
                              use="pairwise.complete.obs", method="pearson"),
            ireg_affect = coef(lm(posaff ~ negaff, na.action=na.exclude))[2])

Look at bivariate IIV metrics data

round(head(istats_bivariate),2)

## # A tibble: 6 × 4
##      id icor_affect icov_affect ireg_affect
##   <dbl>       <dbl>       <dbl>       <dbl>
## 1   203       -0.63       -0.29       -0.98
## 2   204       -0.77       -1.18       -0.7 
## 3   205       -0.44       -0.17       -0.53
## 4   208       -0.35       -0.06       -0.88
## 5   211       -0.08       -0.02       -0.17
## 6   214       -0.85       -0.83       -0.89

Describe the bivariate IIV metrics

describe(istats_bivariate[ ,-1])

##             vars  n  mean   sd median trimmed  mad   min  max range  skew
## icor_affect    1 30 -0.49 0.27  -0.54   -0.50 0.32 -0.90 0.01  0.91  0.30
## icov_affect    2 30 -0.39 0.41  -0.25   -0.32 0.31 -1.73 0.01  1.74 -1.61
## ireg_affect    3 30 -0.58 0.35  -0.62   -0.58 0.43 -1.21 0.03  1.24  0.03
##             kurtosis   se
## icor_affect    -1.24 0.05
## icov_affect     2.28 0.07
## ireg_affect    -1.20 0.06

Correlations among the bivariate IIV metrics

round(cor(istats_bivariate[ ,-1], use="pairwise.complete.obs"),2)

##             icor_affect icov_affect ireg_affect
## icor_affect        1.00        0.73        0.84
## icov_affect        0.73        1.00        0.51
## ireg_affect        0.84        0.51        1.00

Note the overlap between the iCor, iCov, and iReg metrics. This will be relevant when extending the iReg into the multilevel modeling framework - to examine between-person differences in intraindividual covariation.

Merge the bivariate IIV metrics with other person-level data

AMIB_bivariate <- istats_bivariate %>%
  left_join(AMIB_persons, by="id")

Examine relation between Affect Bipolarity and neuroticism

Look at association between bfi_n and icor_affect:

AMIB_bivariate %>%
  ggplot(aes(x=bfi_n, y=icor_affect)) +
  geom_point(colour="gray40") +
  geom_smooth(aes(group=1), method=lm, se=TRUE, fullrange=TRUE, lty=1, 
              linewidth=2, color="#9A23DE") +
  xlab("Neuroticism") + ylab("Affect Bipolarity") +
  theme_classic() +
  theme(axis.title=element_text(size=16),
        axis.text=element_text(size=12),
        plot.title=element_text(size=16, hjust=.5, face="bold")) +
  ggtitle("Between-Person Association Bivariate IIV")

## `geom_smooth()` using formula = 'y ~ x'

Fit model with the neuroticism predicting the individual slope score for the association between positive and negative affect.

reg_bivariate <- lm(icor_affect ~ bfi_n, 
                     data=AMIB_bivariate,
                     na.action = na.exclude) 
summary(reg_bivariate)

## 
## Call:
## lm(formula = icor_affect ~ bfi_n, data = AMIB_bivariate, na.action = na.exclude)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.42499 -0.20524 -0.05017  0.24066  0.51684 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) -0.42619    0.16959  -2.513    0.018 *
## bfi_n       -0.01979    0.04978  -0.398    0.694  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2788 on 28 degrees of freedom
## Multiple R-squared:  0.005616,   Adjusted R-squared:  -0.0299 
## F-statistic: 0.1581 on 1 and 28 DF,  p-value: 0.6939

• For the prototypical person, affect bipolarity score for the sample is -0.43
• For every unit increase in neuroticism, affect bipolarity is expected to change by -0.02 (however, this association is not statistically significantly different from zero p = 0.69). The interpretation is that individuals who have scores higher on neuroticism tend to have stronger affect bipolarity (scores further from zero).

Calculating multivariate IIV metrics

First, set the variables of interest (be sure the id variable is included in the first position)

vars <- c("id", "enthusiastic", "happy", "proud", "relaxed", "calm",
          "sluggish", "sad", "disappointed", "angry", "nervous") 

Calculate the metrics. Each one takes a few steps

# remove ID 218 (due to zero variance in several items)
AMIB_daily_sub <- AMIB_daily[AMIB_daily$id!=218, ]

#calculating the determinant for each individual
istat_idet <- ddply(AMIB_daily_sub[ ,vars], "id",
                  function(x) {
                    det(cor(x[,-1],use="pairwise.complete.obs",
                            method="pearson"))} )
#renaming columns
colnames(istat_idet) <- c("id","idet")

#calculating Frobenius norm for each individual
istat_iFnorm <- ddply(AMIB_daily_sub[ ,vars], "id",
                  function(x) {
                    norm(cor(x[,-1],use="pairwise.complete.obs",
                             method="pearson"),type="F")} )
#renaming columns
colnames(istat_iFnorm) <- c("id","iFnorm")

#merging all together
istats_multivar <- istat_idet %>%
  left_join(istat_iFnorm, by="id")

Look at multivariate IIV metrics data

round(head(istats_multivar),2)

##    id idet iFnorm
## 1 203 0.00   5.64
## 2 204 0.00   5.11
## 3 205 0.03   4.01
## 4 208 0.00   4.43
## 5 211 0.00   4.72
## 6 214 0.00   4.92

Describe the multivariate IIV metrics

describe(istats_multivar[ ,-1])

##        vars  n mean   sd median trimmed  mad  min  max range skew kurtosis   se
## idet      1 29 0.01 0.01   0.00    0.00 0.00 0.00 0.04  0.04 2.00     2.94 0.00
## iFnorm    2 29 4.89 0.65   4.89    4.83 0.56 3.97 6.98  3.01 1.04     1.66 0.12

Correlations among the multivariate IIV metrics

round(cor(istats_multivar[ ,-1], use="pairwise.complete.obs"),2)

##         idet iFnorm
## idet    1.00  -0.64
## iFnorm -0.64   1.00

Merge the multivariate IIV metrics with other person-level data

AMIB_multivariate <- istats_multivar %>%
  left_join(AMIB_persons, by="id")

Examine relation between emotion system variability and neuroticism

Look at association between idet and bfi_n:

AMIB_multivariate %>%
  ggplot(aes(x=bfi_n, y=idet)) +
  geom_point(colour="gray40") +
  geom_smooth(aes(group=1), method="lm", se=TRUE, fullrange=TRUE, lty=1,
              linewidth=2, color="#9A23DE") +
  xlab("Neuroticism") + ylab("Emotion System\nVariability") +
  theme_classic() +
  theme(axis.title=element_text(size=16),
        axis.text=element_text(size=12),
        plot.title=element_text(size=16, hjust=.5)) +
  ggtitle("Between-Person Association Multivariate IIV")

## `geom_smooth()` using formula = 'y ~ x'

And the formal model test using (cross-sectional) regression.

reg_multivariate <- lm(idet ~ bfi_n, 
                     data=AMIB_multivariate,
                     na.action = na.exclude) 
summary(reg_multivariate)

## 
## Call:
## lm(formula = idet ~ bfi_n, data = AMIB_multivariate, na.action = na.exclude)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -0.0064492 -0.0046937 -0.0038566  0.0003153  0.0305352 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)
## (Intercept)  0.0082545  0.0064399   1.282    0.211
## bfi_n       -0.0009027  0.0018611  -0.485    0.632
## 
## Residual standard error: 0.009515 on 27 degrees of freedom
## Multiple R-squared:  0.008638,   Adjusted R-squared:  -0.02808 
## F-statistic: 0.2353 on 1 and 27 DF,  p-value: 0.6316

• For the prototypical person, emotion system variability score for the sample is 0.01
• Individuals who are higher on neuroticism tend to have lower emotion system variability. For every unit increase in neuroticism, emotion system variability is expected to change by 0 (however, this association is not statistically significantly different from zero p = 0.63).

Conclusion

This tutorial illustrated a general approach for calculting intraindividual variability metrics for univariate, bivariate, and multivariate data - metrics that can be used as indicators of specific constructs (e.g., lability, complexity, etc.).

We hope this foundational material prompts consideration of a larger set of intraindivdual variability constructs and the metrics used to operationalize those constructs.

Have fun exploring!