Using Lavaan To Calculate Alpha And Omega

[stephenba...@gmail.com]

Hi,

I am trying to run omega and alpha statistics using a combination of the lavaan package and the SEMTools package
I am a complete newbie on this so apologies if the question seems basic. I was inspired by the post https://drsimonj.svbtle.com/how-to-calculate-internal-consistency
However a lot of the litrature seems to specify that boot strapping and confidence intervals should be used when calculating Cronbach Or Omega

I have 3 constructs taken from the `lavaan` data HolzingerSwineford1939
My aim is to run both Alpha and Omega on my constructs that i have created with the final result being to plot the point statistic with error bars showing the confidence interval of each construct.

So in theory i want get back my alpha/omega statistic and the confidence interval so i can plot my point with error bars

I am getting the below error when i run it

Warning messages:
1: In lav_object_post_check(lavobject) :
  lavaan WARNING: some estimated ov variances are negative
2: In lav_object_post_check(lavobject) :
  lavaan WARNING: some estimated ov variances are negative
3: In lav_object_post_check(lavobject) :
  lavaan WARNING: some estimated ov variances are negative
4: In lav_object_post_check(lavobject) :
  lavaan WARNING: some estimated ov variances are negative

Can anyone help me how to trouble shoot?
Below is the code which i appreciate needs a lot of work :)

#------------------------
# CODE
#-------------------------
library(lavann)
library(tidyverse)
library(semTools)

HS.model <- ' visual  =~ x1 + x2 + x3
              textual =~ x4 + x5 + x6
              speed   =~ x7 + x8 + x9 '


# Fit the model
boot.fn <- function(model,data, indices, testtype) {
 
  library(semTools)
  d <- data[indices,]
 
  if(testtype == 'alpha'){
    slice_index = 1
  } else {
    slice_index = 2
  }
 
  fit <- cfa(model, d)
  # Slicing either 1 or 2 will result in either Omega or Alpha results being returned
  rel <- as.data.frame(reliability(fit)) %>%
    slice(slice_index)
 
  rel <- as.numeric(rel)
 
  return(rel)
 
}

results.alpha <- boot(data=HolzingerSwineford1939, statistic=boot.fn, R=1000, model = HS.model, testtype = 'alpha')
results.omega <- boot(data=HolzingerSwineford1939, statistic=boot.fn, R=1000, model = HS.model, testtype = 'omega')

boot.ci(results.alpha)


Im taking my information on how to create boot strapping from 'R In Action' Pg 306 and ISLR page 196

Cheers
Steve

[Terrence Jorgensen]

I am getting the below error when i run it. Can anyone help me how to trouble shoot?

Those aren't errors, they are warnings.  When you have a variance estimate close to zero, sampling error across bootstrap samples will almost certainly result in estimates in some samples being negative.  Did you get that warning when you fit your model to the original data?  

If you still get results, I wouldn't worry about these warnings.  The purpose of bootstrapping is to simulate a sampling distribution, not for every single estimated variance to be within theoretical bounds. 

Below is the code which i appreciate needs a lot of work :)

Looks good to me.  Two notes on possibilities:

1. Bootstrapping takes time.  You can write your custom function to return a vector (both alpha and omega) instead of only one number, that way you only have to bootstrap once.  The boot.ci() function has an argument  to select the index of the parameter you want a CI for.  
2. You should also select the "bca" type, which tends to provide better coverage in non-normal sampling distributions.

Terrence D. Jorgensen

Postdoctoral Researcher, Methods and Statistics

Research Institute for Child Development and Education, the University of Amsterdam

UvA web page: http://www.uva.nl/profile/t.d.jorgensen

[stephenba...@gmail.com]

Hi Terence

Thank you very much for your help. This problem has been gnawing at me for a week :)
I have updated my code as you recommended and its perfect although it seems i have to put one construct at a time through it which is ok of course.
I have learned a lot :)

# Fit the model
boot.fn <- function(model,data, indices) {

 
  library(semTools)
  d <- data[indices,]
 

  fit <- cfa(model, d)
  rel <- as.data.frame(reliability(fit)) %>%      
    slice(1:2) %>%
    as.vector()
 
   rel <- as.numeric(return(t(rel)))
 
  return(t(rel))
 
}
results<- boot(data=HolzingerSwineford1939, statistic=boot.fn, R=1000, model = HS.model)

alpha.ci <- boot.ci(results,type = 'bca', index = 1)
omega.ci <- boot.ci(results,type = 'bca', index = 2)

I should now be able to make a dataframe with confidence intervals of the scores for each construct

Thanks again for your help
Cheers
Steve