Heteroscedastic one-factor models for Visual Analogue Scales

(1)

*Correspondence should be addressed to Esther A. Lietaert Peerbolte, Pieterskerkhof 11, 2311 SP, Leiden, The Netherlands (e-mail: estherlietaertpeerbolte@gmail.com).

Heteroscedastic one-factor models for Visual Analogue Scales

Masterthesis by Esther A. Lietaert Peerbolte 1_*

Supervisors: Dylan Molenaar 1_{& Gideon J. Mellenbergh} 1

1_{Department of Psychology, Faculty of Social and Behavioral Sciences, University of Amsterdam,} Weesperplein 4, 1018 XA Amsterdam, the Netherlands

Abstract

Homoscedasticity is an implicit assumption in factor analysis using maximum likelihood estimation. It implies that the error variance and common factor are unrelated. Violations, known as heteroscedasticity, are difficult to detect and often exist when Visual Analogue Scales are used. While previous studies have introduced heteroscedastic one-factor models they only focused on a monotonically increasing or decreasing relation between the error variance and common factor, a relation that is often not appropriate for responses on Visual Analogue Scales. Therefore this study focuses on finding a better method to fit a one-factor model on this type of data, using both data-based and model-data-based methods. First, the arcsine transformation was considered with a standard one-factor model. This transformation did not provide better results than the fit of the factor model on the raw data and is therefore rejected. Second, the heteroscedastic model presented in Hessen and Dolan (2007) and Molenaar et al (2010) was extended using the beta function to relate the error variance to the common factor. Simulations showed that this model cannot handle severely non-normal data as is often the case with Visual Analogue Scales. The model was therefore rejected. Finally, a mixture model was created in which all items measure one latent trait and responses are influenced by error and the latent trait which both follow the beta function. Simulations and an application on real data showed the model is better suited for usage on data measured using line segments than the standard one-factor model. In conclusion, the mixture model presented in this study is a suitable method to fit a one-factor model to data measured on line segments and is therefore preferred over the standard one-factor model.

(2)

1. Introduction

Factor analysis is a latent trait model for continuous items. Using linear regression it searches for the number of latent variables that best explain the correlation structure of the manifest variables. By explaining these correlations the factor structure of the data is revealed (Lawley & Maxwell, 1971). The discovery of such structures has lead to valuable information, for instance about the structure of personality (e.g. Digman, 1990; Goldberg, 1990) and the structures of psychological disorders (e.g. Nelson et al, 1999). As for any other statistical method, linear factor analysis has its own sets of assumptions: (1) item responses are explained by a number of variables, of which at least one is latent, (2) local independence between the items, (3) a linear relation between the manifest and latent variables, and (4) homoscedasticity of the error variance (Mellenbergh, 1994).

While for the first three assumptions violations can be detected and accommodated in the factor model, the assumption of homoscedasticity is difficult to detect and often violated when Likert scale items or line segments are used, possibly resulting in a lack of model fit (Mellenbergh, 1994). Homoscedasticity is the assumption of homogeneous variances. It means that the precision of measurement, or error variance, and the values of the latent variable are independent of each other. In other words, the variance of the measurement error is approximately the same at every level of the latent variable. When the assumption is violated and there is any sort of relation between the latent variable and the residual variance, this is called heteroscedasticity. Figure 1 illustrates the difference between homoscedasticity and heteroscedasticity. In Figure 1a it is clearly visible that the amount of scatter at each level of the latent variable is approximately the same, indicating homoscedasticity. However, in Figure 1b the amount of scatter becomes larger as the latent variable increases. This implies linear heteroscedasticity.

This assumption of homoscedasticity is the result of a frequently used method to find the best latent structure that fits the data: maximum likelihood (ML, Lawley, 1943). This method requires the observed scores to be distributed according to a multivariate normal distribution. This assumption implies that (1) the latent variable is normally distributed, (2) the factor loadings are equal for all respondents, and (3) the error is homoscedastic and normally distributed. Thus, the fact that the residuals are often heteroscedastic when maximum likelihood is used is a problem that requires a solution: heteroscedastic models.

However, even though the assumption of homoscedasticity is often violated only a few studies have focused on heteroscedastic factor models (i.e. Bollen, 1996; Lewin-Koh and Amemiya, 2003; Meijer and Mooijaart, 1996). Hessen and Dolan (2007) created a heteroscedastic one-factor model to create the ability to test for homoscedasticity and the possibility to adequately fit models that violate the assumption of homoscedasticity. They did so using marginal maximum likelihood

(3)

estimation. By turning the residual variance into an exponential function and thereby relating it to the latent variable they were able to create a heteroscedastic one-factor model.

As mentioned earlier, this assumption of multivariate normality, however, not only implies homoscedasticity, but also implies that the latent variable is normally distributed and that the factor loadings are equal for all individuals. Molenaar, Dolan and Verhelst (2010) extended the previous heteroscedastic model by also dropping these two assumptions. Their analyses show that a non-normal latent distribution can be estimated together with heteroscedasticity and that factor loadings dependent on the common factor and heteroscedasticity can also be estimated together, but that non-normality of the latent variable and level dependent factor loadings cannot be estimated together. This model provides explicit tests of these assumptions instead of an omnibus test for all three assumptions combined and the ability to use marginal maximum likelihood estimation in a one-factor model when its assumptions are violated.

(4)

However, these previous heteroscedastic models are limited in two ways. First, these models do not take into account data with a bounded interval. This type of violation is easily accommodated as most datasets can be transformed to fit this interval. Second, the relation between the error variance and the latent trait was linear in these models. The variance for example increased as the latent trait increased, as was shown in Figure 1b. While this might work for some data, it will not be the appropriate model for all sorts of data. It is therefore necessary to create heteroscedastic models that are more flexible. For line segments it is for instance possible that error variances decrease near the extremes of the scale (Mellenbergh, 1994). For these reasons this study will focus on possible one-factor models for data that can be severely non-normal. This would be especially beneficial for data measured using line segments.

Line segments, also known as Visual Analogue Scales (VAS), are a popular method to assess for example pain (e.g. Hawker et al, 2005; Hjermstad, 2011; Price et al, 1983). In a closed response situation as shown in Figure 2 respondents are allowed to check any point on the line in order to indicate the severity of their perceived pain (Samijima, 1973). The scale in Figure 2 is 9 centimeters long, resulting in a range from 0 mm (no sensation) to 90 mm (the most intense sensation imaginable). The patient is asked to mark the point on the line that corresponds to his or her experienced severity of pain. In this case the respondent has answered 3 centimeters to the right of “no sensation”, resulting in a score of 30.

Figure 2. An example of a response on the real line, otherwise called a visual analogue scale (VAS).

These line segments have two points of objective certainty: the boundaries, since they are labeled and represent the extremes of the item tested. However, everything in between is sensitive to subjectivity. A score of 30/90 for one person may mean something different for someone else and is therefore to a certain extent subjective. The closer you get to the labeled boundaries the more certain the interpretation of those parts of the line become and therefore the more objective they become. At the same time the boundaries also mark the end of the line with the consequence that not much variance can exist near those points. Thus, it is likely that less error variance will be found

(5)

around the boundaries, whilst more error variance will be found in the middle of the line segment, taking for example the shape of a parabola.

In order to model the error variance in such a way that it resembles a parabola Hessen and Dolan (2007) proposed to use a polynomial. However, a polynomial is very limiting. A polynomial indicates that the largest or smallest residual variance is always exactly in the middle of the scale. In addition, it indicates that the function between the common factor and the variance is always symmetrical, which is not necessarily realistic for responses on line segments. Thus, these assumptions are too strict and a more appropriate function needs to be used.

For these reasons this study will focus on possible heteroscedastic one-factor models that can model severe non-normality. Three methods are explored in this paper. First, a data-based approach will be tested. This will be done because it would be a much easier and therefore perhaps more often used method for line segment data than a new statistical model. In the second section the arcsine transformation, a variance stabilizing transformation, will therefore be used to see if it provides a better fit than the raw data when fitting a standard one-factor model. If this does not provide a better fit, two different models will be created and tested using simulations and real data. Second, the beta variance model will be tested in the third section. It will have a normally distributed latent variable with error variances that follow the beta function. This resembles closely the models of Hessen and Dolan (2007) and Molenaar, Dolan and Verhelst (2010) with the exception that the logarithm of the error variance is no longer linearly related to the common factor. This model is more appropriate than the suggested polynomial as it is more flexible. However, this model is limited in that it does not handle extremely skewed data very well. Finally, a second model will be tested that is better equipped to deal with such data. This is because the model is a mixture model with a latent variable that follows the beta function, in which the shape parameters are equal to form a symmetric distribution, and an item-specific component that also follows a beta distribution, but is free in its parameter values, which will model the heteroscedasticity. This model will be discussed in the fourth section. Finally, a general discussion will follow in the fifth section.

2. Arcsine Transformation

As mentioned earlier a data-based solution to the problem of heteroscedasticity in line segment data might be more beneficial than a model-based solution as the former is easier to use and understand. Therefore before any models are created the option of transforming data will be explored. The arcsine transformation is chosen as the most appropriate transformation for three reasons: (1) it makes skewed data more normally distributed, (2) it can only be used on data on the [0,1] interval, which, as explained in the introduction, is the chosen interval for this study and (3) it

(6)

makes the variance more homogeneous (Novick and Jackson, 1974). The transformation is formulated as

𝑔 = sin−1√𝑥, (1)

where g is the transformed data and x is the observed data. Because of this function the interval for the transformed data is [0, π/2], where the maximum is approximately 1.57.

2.1. Application

In order to test the efficiency of the arcsine transformation it will be applied to a dataset. As the main goal of this article is to find a better method to fit a one-factor model to data measured using Visual Analogue Scales the only requirement for an appropriate dataset is that responses are measured on line segments. In this case a version of the Bermond-Vorst Alexithymia Questionnaire (BVAQ; Vorst & Bermond, 2011) was used.

2.1.1. Background

Alexithymia stems from the Greek words ἀλέξις θυμός, translating to “no words for emotions”. It is a personality trait characterized by difficulty in identifying and describing emotions (Taylor, 1984). Individuals with this trait think in a concrete and reality-based manner that is more externally oriented. Because of their restricted imagination they have few dreams and fantasies. While alexithymia is officially classified as a personality trait and not as a mental disorder, it is linked with an increased risk of developing medical or psychiatric disorders. At the same time it also appears to reduce the likelihood of responding to the conventional treatments of these disorders (Haviland, Warren & Riggs, 2000). For these reasons alexithymia has been of interest in both personality and clinical research.

The BVAQ consists of 40 items measuring an affective and a cognitive component of alexithymia. The affective component consists of two subscales, emotionalizing and fantasizing, whilst the cognitive component has three subscales, verbalizing, identifying and analyzing. Each of these subscales are measured with 8 items. While these items are usually scored on a 5-point scale another version of the questionnaire was used for this project, in which responses are measured on a line segment ranging from 0 to 100. Since this study focuses on scales with a [0,1] interval, all responses were divided by 100. Half of the items were contra-indicative and therefore recoded before analyses were conducted. A high score on an item or scale indicates a high tendency towards alexithymia. Because this project focuses solely on a one-factor model each subscale will be tested separately.

(7)

2.1.2. Results

The line segment-version of the BVAQ was completed by 413 first year psychology students. The reliability was calculated for each of the subscales. Cronbach’s alpha values are 0.794, 0.859, 0.857, 0.791 and 0.824 for emotionalizing, fantasizing, verbalizing, indentifying and analyzing, respectively, indicating high reliability for each of the subscales. Therefore all five subscales are deemed appropriate for use in this application study. Means, standard deviations and ranges for all five scales are presented in Table 1.

A baseline was created by fitting the standard Maximum Likelihood one-factor model to each of the subscales before a transformation was applied. The Root Mean Square Error Approximation (RMSEA) and Root Mean Square Residual (RMSR) are presented in Table 2 as absolute goodness of fit measures. Based on these results the fit of the standard one-factor for the five subscales ranges from poor to mediocre, with Fantasizing arguably being the only convincingly acceptable fit.

Table 1. Means, standard deviations and ranges of the variables in the application

Verbalizing Emotionalizing Fantasizing Identifying Analyzing

Mean .50 .41 .49 .38 .37

SD .27 .25 .29 .24 .24

Min 0 0 0 0 0

Max 1 1 1 1 1

Table 2. Goodness of Fit measures of the standard one-factor model on the raw data for all five

subscales

RMSEA .132 .128 .068 .119 .088

RMSR .069 .085 .039 .088 .054

Next, the data were transformed using Equation 1. The same one-factor model that was applied to the baseline was now applied to the transformed data. RMSEA and RMSR values are presented in Table 3. Fantasizing has arguably the best fit, but overall the fit is still mediocre to poor. When comparing Table 2 and 3 it is clearly visible that the differences between the goodness of fit measures are minimal. In most cases they are a little better when the raw data is used and sometimes they are unchanged. Only the RMSR of Fantasizing improved by transforming the data and this is a difference of 0.003. Therefore it can be concluded that using the arcsine transformation on the data does not provide a better fit of the standard Maximum Likelihood one-factor model on line segment data compared to the fit of the same model on the observed data.

(8)

Table 3. Goodness of Fit measured of the standard one-factor model on the transformed data of all

five subscales

RMSEA .137 .135 .068 .138 .100

RMSR .071 .086 .036 .101 .058

2.2. Conclusion

The arcsine transformation is thus dismissed as a solution to the problem of heteroscedasticity in line segment data. Since the arcsine transformation does not even seem promising to use it, it seems better to test a model-based solution instead. Therefore the rest of the study will focus on the creation of heteroscedastic one-factor models.

3. Beta variance model

The one-factor model for person i and item j is formulated as

𝑦𝑖𝑗 = 𝜐𝑗+ 𝜆𝑗𝜂𝑖+ 𝜀𝑖𝑗, (2)

where yij is the observed score, υj is the intercept, λj is the factor loading relating the item to the common factor, ηi, and εij is the residual or error parameter. In Hessen and Dolan (2007) and Molenaar et al (2010) the residual variance was modeled by using the simple regression equation

log(𝜎2𝜀𝑗|𝜂) = 𝛽𝑗0+ 𝛽𝑗1𝜂𝑖 , (3)

which is referred to as the minimal heteroscedastic model. In this model 𝛽𝑗1 is the heteroscedasticity parameter. The parameter 𝛽𝑗0 accounts for the residual variance in the observed variable, independent of the common factor.

As mentioned earlier the main difference between this study and Hessen and Dolan (2007) and Molenaar et al (2010) is that in these studies the logarithm of the residual variance was modeled to be linearly related to the latent factor, taking into account skewness, while this model should be better equipped to deal with severe non-normality of the data. While Hessen and Dolan (2007) proposed a polynomial to do so, its implicit assumptions are too strict for data measured using line segments.

(9)

The beta function solves this problem and therefore seems better suited to model the residual variance in Visual Analogue Scales, as it allows not just for skewness, but also kurtosis. The beta function without the normalizing constant is:

𝑥𝑖𝛼−1(1 − 𝑥𝑖)𝛽−1, (4)

where xi represents the observed score for person I and α and β are the shape parameters. This function is defined on interval [0,1] and has two positive shape parameters, α and β. When the two shape parameters are equal the shape of the function resembles a Gaussian distribution. Using the beta function gives the advantage that responses are bounded between 0 and 1 and that the shape parameters can be modeled to create functions of different shapes, such as bimodal functions (0 < α < 1, 0 < β < 1), uniform functions (α=β=1), symmetrical functions (α = β, α > 1, β > 1), and skewed functions (α ≠β, α > 1, β > 1). Examples of such functions can be found in Figure 3.

Figure 3. Examples of the beta function with different sets of shape parameters.

The function defined in Equation 4 is therefore used to relate the variance and the common factor in this model. A baseline parameter was added to the function as a way to model the amount of variance in the case of homoscedasticity (α=β=1). This creates the conditional residual variance function

(10)

𝜎2𝜀𝑗 | 𝜂 = 𝛽𝑗0 ( 𝜂𝑖𝛼𝑗−1 (1 − 𝜂𝑖)𝛽𝑗−1) , (5)

where β0j is the baseline parameter.

As this model can be seen as an extension of the previously mentioned models of heteroscedasticity the current model is the same as the one in Molenaar et al (2010) with the exception that the conditional residual variance function is changed from Equation 2 to Equation 5 and the mean and standard deviation of the common factor have different values. The reason why the values are different is because data in this model will always have a minimum of 0 and a maximum of 1. Therefore a mean of 0 and standard deviation of 1 are not the best choice. This is why they are changed to 0.5 and 0.1, respectively. Further details of the model can be found in Molenaar et al (2010).

3.1. Simulation

The viability of the model is tested with a simulation study. In this study heteroscedastic data will be created using the parameters described above. The values of these parameters will then be estimated under the model.

3.1.1. Design

I simulated data for a one-factor model with heteroscedasticity using Equations 2 and 5. The shape parameters of the residual variance were chosen based on real line segment data. The dataset is described in Section 4.3 and plots are presented in Appendix D. Nearly all of the items in the dataset were either slightly or extremely skewed. Some resembled a normal distribution and a few a bimodal distribution. Therefore the shape parameters were manipulated in such a way that all of these conditions were represented (4 levels of heteroscedasticity: [α,β] = [0.9, 0.9], [3,4], [4,4], [7,4]; 1 level of homoscedasticity: [1,1]). The number of simulees was manipulated as well (3 levels: 300, 500, 700 simulees), creating 15 conditions (5 x 3). In order to create data between 0 and 1, a baseline parameter value, βj0, was chosen for each of the five levels of residual variance. Because the baseline parameter was often too small for the model to estimate its logarithm was used instead. For the same reason as the baseline parameter was introduced the intercept, υi, and factor loading, λi, were chosen for every fifth item.A total of 15 items were simulated, resulting in the same intercept and factor loading for three items. The corresponding levels of βj0 are shown in Table 4. The intercepts and factor loadings for the items are shown in Table 5. The plots of the five levels of residual variance are shown in Figure 4. The factor scores were randomly drawn from a normal distribution with mean 0.5 and standard deviation 0.1.

(11)

Table 4. Selected shape parameter conditions and their corresponding height parameters. α β βj0 log(βj0) 0.9 0.9 0.0009 -7.01 3 4 0.1 -2.30 4 4 0.4 -0.92 7 4 0.1 -2.30 1 1 0.0009 -7.01

Table 5. Selected parameter values for the intercept and factor loadings for every condition.

Item i υj λj 1, 6, 11 0.04 0.9 2, 7, 12 0.02 0.93 3, 8, 13 0.06 0.92 4, 9, 14 0.03 0.94 5, 10, 15 0.05 0.95

Figure 4. Plots for the residual variance.

Several restrictions were imposed for identification purposes. First, the mean and standard deviation of the underlying factor were fixed to 0.5 and 0.1, respectively. Second, because the shape parameters can only be positive the lower bound of both was set to 0.01. Last, the initial values for the scale parameters were set to 2, while the remainder of the initial values for estimation

(12)

parameters were set to 0. The reason why a different initial value was chosen for the scale parameters compared to all other parameters is that it is the first rounded value the parameters can take when they are approximately normally distributed.

For each of the 15 conditions 100 datasets were simulated. The model was then fitted to these datasets using OpenMX (Boker et al, 2011). The code can be found in Appendix A. In all analyses we will use 30 quadrature points in order to keep the model as close as possible to the model in Molenaar et al (2010).

3.1.2. Results

In order to establish the goodness of fit Root Mean Square Error (RMSE) values were calculated. The RMSE is the square root of the variance of the residuals, otherwise known as the standard deviation of the unexplained variance. As it compares the true value of the parameter to the estimated value the closer the RMSE gets to zero the better the parameter recovery. The results of the error variance parameter recovery are shown averaged per condition in Table 6 and the results of the estimation of the other parameters, the intercept and factor loadings of the data, are shown averaged per three items in Table 11 in Appendix D.

Visual inspection of these results show as expected the estimation of all parameters improves as the number of simulees increases. The intercept and factor loadings are estimated rather unbiased, but the estimates of the baseline parameter and shape parameters of the residual variance are more biased. The model appears to handle homoscedastic (α = β = 1) data and data that is almost homoscedastic (α = β = 0.9) better than heteroscedastic data. These two conditions are estimated rather unbiased, while the model could for example not estimate the residual variance parameters of condition [7,4] well with 300 simulees and could not estimate them at all with 500 and 700 simulees. R returns a “non-zero status code 6”, also known as “Mx status RED”, implying that “the model does not satisfy the first-order optimality conditions to the required accuracy, and no improved point for the merit function could be found during the final linesearch”. Results for these conditions are therefore omitted in Table 6. When the data is symmetrically distributed as in condition [4,4] the model seems to work better than when the data is slightly skewed as in condition [3,4].

3.3. Conclusions

The results of the simulation indicate that the model works on homoscedastic data and symmetrically distributed data, but cannot handle heteroscedastic data well. While the model works well enough when the data is slightly skewed, it does not function appropriately when the data is more extreme. The model is thus promising as it is the first heteroscedastic factor model that can

(13)

handle slightly skewed data. However, because this study focuses mostly on severely skewed data the model is dismissed as a solution to the problem of heteroscedasticity in line segment data and a new model will be tested.

The beta function seemed the best possible function to model the residual variance in an extension of the previous heteroscedastic models. Since it did not function as well as expected, it can be concluded that the previous models cannot be properly adapted to work on data measured using Visual Analogue Scales.

(14)

Table 6.

True values, mean observed values and RMSE values of estimated parameters in the residual variance

Βo α β

N (α,β) Log(True) Obs (SD) RMSE True Obs (SD) RMSE True Obs (SD) RMSE

300 (0.9, 0.9) -7.01 -6.40 (.139) 0.063 0.9 1.29 (.096) 0.398 0.9 1.36 (.107) 0.467 (1, 1) -7.01 -6.54 (.121) 0.489 1 1.28 (.083) 0.293 1 1.36 (.087) 0.372 (3, 4) -2.30 -0.59 (.299) 1.733 3 3.93 (.191) 0.947 4 5.72 (.246) 1.737 (4, 4) -.092 -0.63 (.214) 0.644 4 4.20 (.146) 0.440 4 4.21 (.157) 0.490 (7, 4) -2.30 -7.38 (.051) 5.077 7 1.50 (.033) 5.499 4 1.25 (.052) 1.754 500 (0.9, 0.9) -7.01 -6.63 (.091) 0.388 0.9 1.14 (.064) 0.245 0.9 1.15 (.069) 0.262 (1, 1) -7.01 -6.76 (.090) 0.268 1 1.14 (.058) 0.155 1 1.17 (.073) 0.180 (3, 4) -2.30 -1.38 (.228) 0.948 3 3.50 (.153) 0.523 4 4.91 (.168) 0.929 (4, 4) -0.92 -0.63 (.214) 0.352 4 4.20 (.146) 0.244 4 4.21 (.157) 0.257 (7, 4) -2.30 - - 7 - - 4 - - 700 (0.9, 0.9) -7.01 -6.76 (.079) 0.268 0.9 1.06 (.057) 0.168 0.9 1.07 (.053) 0.174 (1, 1) -7.01 -6.81 (.103) 0.229 1 1.11 (.067) 0.131 1 1.13 (.074) 0.146 (3, 4) -2.30 -1.78 (.158) 0.557 3 3.28 (.111) 0.302 4 4.54 (.113) 0.551 (4, 4) -0.92 -0.77 (.178) 0.226 4 4.10 (.123) 0.154 4 4.11 (.125) 0.164 (7, 4) -2.30 - - 7 - - 4 - -

(15)

4. Beta Mixture Model

Because the beta variance model could not handle heteroscedastic data and severe non-normality in the data well enough another model is needed. A suited option seems a mixture model. Mixture models are probabilistic models that take into account subgroups in the overall population. In other words they can handle heterogeneity in the population by identifying homogeneous groups, also known as latent classes. As a consequence they handle extreme data rather well.

This beta mixture model is a Bayesian latent mixture with one latent factor, a common factor that is the same for everyone and for every item, and an item-specific component. As a result the answers on items measuring the same underlying construct can be very differently distributed. The model is therefore very flexible and will be better equipped to deal with extreme data than the last model.

The model can be formulated as

𝑓(𝑦) = ω𝑔(𝜂) + (1 − ω)ℎ(𝜀), (6)

where f(y) is the observed distribution, 𝜔 is the mixing proportion that can range from 0 to 1, the common factor is

𝜂 ~ 𝛣(𝛼, 𝛼), (7)

where B(𝛼,𝛼) is a Beta function with shape parameters that are equal in order to make the distribution symmetrical so it is comparable to other models in which the latent trait is defined as a normal distribution, and the error term is

𝜀 ~ 𝛣(𝛽, 𝛾). (8)

As is shown here, the function of the mixture model can be cut in two halves; the first half being the part of the common factor and the second half being the part of the error variance. The common factor and the item-specific component are related to each other through the mixing proportion. When the mixing proportion is small the first part of the equation is very small, while the second part of the equation is large. This essentially means that a person’s score on an item is based more on that item than on their latent score of the construct that is being measured. Therefore the second part of the equation can be seen as an error term. In contrast when the mixing proportion is large the person’s score is based mostly on their latent score, indicating that the item is rather

(16)

unbiased. The mixing proportion is thus relatable to the factor loadings in the beta variance model. It is also relevant for item analysis as a low mixing proportion indicates a biased item.

4.1. Simulation

Before we will apply the model to real data we will first test its viability in a simulation study.

4.1.1. Design

I simulated data for a one-factor model with heteroscedasticity using Equation 6. As mentioned earlier the first half of the equation represents the latent factor and the second half represents the error term. The common factor was simulated by using the beta function with shape parameters which were always equal to each other in order to create a symmetrically distributed latent factor, but varied over the conditions. The error values were also created using the beta function but the shape parameter values created different sorts of distributions and thus varied over the conditions. The mixing proportion was always set equally for all items, but varied over the conditions as well. All these parameter values were chosen based on the same dataset as the beta variance model and are presented in Table 7. Plots for each of the combinations are shown in Figure 5. Five different conditions were created by combining these parameters. Every condition was simulated for 300, 500 and 700 simulees, creating a total of 15 (5x3) conditions.

Table 7. Selected shape parameter conditions for both the latent variable and the error term, and

their corresponding mixing proportion values.

ω α β γ 0.2 3 1 10 0.6 3 2 7 0.5 2 5 2 0.7 3 7 1 0.5 2 9 2

As mentioned earlier a restriction that is posed in this model is that the shape parameters of the latent factor are equal to each other in order to create a symmetrically shaped distribution. Another restriction set to the model is that the shape parameters of both the latent factor and the error must always be at least 0.01 in order to ensure they are positive values. The initial values for the shape parameters were set to 2 once again because it is the first rounded value the parameters can take when they are approximately normally distributed. The initial value of the mixing proportion is set to 0.5, because it means the common factor and error are equally weighted.

(17)

Figure 5. Plots for the different shape parameter conditions as presented in Table 8.

For each of the 15 conditions 100 datasets were simulated. Values of 1 were changed to 0.99 and values of 0 were changed to 0.01 as the model was not able to handle the values of the boundaries. The model was fitted to these datasets using OpenBUGS (Thomas et al, 2006). For each simulation a total of 10.000 samples were drawn of which half were burned. This was done to avoid local minima and maxima. The R-code and OpenBUGS code can be found in Appendix B and C, respectively.

4.1.2. Results

As for the beta variance model RMSE values were computed in order to establish the parameter recovery of the model. These values combined with the true and mean estimated values are presented in Table 8. What can be seen is that, as expected, the estimation becomes more unbiased as the number of simulees increases. What is also apparent is that the mixing proportion and shape parameter of the latent variable are estimated well in all conditions. As the true values of the shape parameters of the items have a much larger range than the previously mentioned parameters, more specifically from 1 to 12, the accuracy of the parameter recovery varies. Small values are estimated rather unbiased, while larger values result in large standard deviations and therefore larger RMSE values. However, these estimates are still acceptable.

(18)

Overall, the results of the simulation for the beta mixture model are satisfying. Values of the mixing proportion and smaller shape parameters are estimated well and show very good goodness of fit values. On top of that, the model never failed to function, like the beta variance model did in the highly skewed condition. Even though the parameter recovery of more extreme values was not perfect, they were still estimated acceptably well. Therefore the mixture model is seen as a suitable one-factor model for heteroscedastic data measured using line segments.

(19)

Table 8.

True values, mean observed values and RMSE values of estimated parameters

ω α β γ

N (β,γ) True Obs (SD) RMSE True Obs (SD) RMSE True Obs (SD) RMSE True Obs (SD) RMSE

300 (1, 10) 0.2 0.21 (0.293) 0.031 3 3.07 (0.336) 0.341 1 1.20 (0.085) 0.213 10 12.35 (1.64) 2.867 (2, 7) 0.6 0.47 (0.098) 0.159 3 2.84 (0.355) 0.389 2 1.91 (0.331) 0.342 7 5.91 (2.197) 2.449 (5, 2) 0.5 0.47 (0.090) 0.093 2 2.11 (0.108) 0.156 5 5.40 (1.427) 1.478 2 2.05 (0.326) 0.330 (7, 1) 0.7 0.65 (0.069) 0.085 3 3.11 (0.329) 0.344 7 7.90 (3.116) 3.233 1 1.12 (0.154) 0.197 (9, 2) 0.5 0.48 (0.055) 0.057 2 2.11 (0.108) 0.153 9 9.16 (2.009) 2.013 2 2.02 (0.282) 0.282 500 (1, 10) 0.2 0.21 (0.023) 0.025 3 3.07 (0.262) 0.271 1 1.20 (0.066) 0.211 10 12.46 (1.238) 2.750 (2, 7) 0.6 0.50 (0.098) 0.140 3 2.86 (0.259) 0.293 2 1.92 (0.274) 0.284 7 6.10 (1.931) 2.130 (5, 2) 0.5 0.49 (0.075) 0.075 2 2.10 (0.113) 0.151 5 5.45 (1.238) 1.315 2 2.06 (0.277) 0.283 (7, 1) 0.7 0.69 (0.055) 0.056 3 3.02 (0.217) 0.216 7 8.91 (2.847) 3.419 1 1.17 (0.135) 0.214 (9, 2) 0.5 0.49 (0.041) 0.041 2 2.07 (0.072) 0.102 9 9.18 (1.544) 1.559 2 2.02 (0.218) 0.220 700 (1, 10) 0.2 0.21 (0.019) 0.023 3 3.02 (0.224) 0.224 1 1.20 (0.057) 0.212 10 12.53 (1.096) 2.756 (2, 7) 0.6 0.51 (0.103) 0.134 3 2.91 (0.231) 0.247 2 1.92 (0.240) 0.252 7 6.16 (1.827) 2.013 (5, 2) 0.5 0.51 (0.068) 0.068 2 2.10 (0.097) 0.137 5 5.44 (1.100) 1.182 2 2.04 (0.226) 0.230 (7, 1) 0.7 0.70 (0.043) 0.043 3 3.00 (0.196) 0.195 7 9.09 (2.411) 3.184 1 1.18 (0.116) 0.211 (9, 2) 0.5 0.50 (0.035) 0.035 2 2.07 (0.078) 0.102 9 9.21 (1.283) 1.300 2 2.02 (0.182) 0.184

(20)

4.2. Application

The model has proven to be successful in the simulation study. Therefore it will now be applied to real data in order to see if it can handle real data sets also. A data set containing items of the Adjective Checklist (ACL) will be used. The alexithymia dataset mentioned in Section 2 will then be used as an extra check.

4.2.1. Background

The ACL measures the presence or absence of 300 adjectives on psychological traits. These traits combined form 30 scales. They are often used to measure views of the self, creativity or stereotypes. In its most common form respondents check the traits that best describe them or someone else. However, in another form all adjectives are measured as line segments. In this case the latter form was used with 6 centimeter long lines. Therefore any adjective is scored from 0 to 60 mm. Examples of items are “aggressive”, “confused”, and “stabile”.

4.2.2. Results

For this study 244 psychology freshmen filled in 22 of the 30 scales with 10 items each, creating a total of 220 items. Of these 22 scales the 5 with the highest reliability were chosen for the application study. These scales are intraception, ideal self, affiliation, order and deference with a reliability of 0.827, 0.803, 0.802, 0.776, 0.738 and, respectively. These Cronbach’s alpha values indicates that all five scales have high reliability. Values are divided by 60 to create a [0,1] interval. As was explained earlier values of 1 are changed to 0.99 and values of 0 are changed to 0.01. Means, standard deviations and ranges for all five scales are presented in Table 9. All estimated parameters of the five scales can be found in Table 10.

Table 11 shows the estimated parameters for the two items with the largest and the two items with the smallest mixing proportion. These represent, respectively, the most unbiased and biased items of all five scales. Figure 6 shows the representing plots. The most biased items are shown in red and the most unbiased items in green. Interestingly enough the estimates and plots of the biased items do not seem significantly different from those of the unbiased items.

Table 9. Means, standard deviations and ranges of the variables in the application study.

Intraception Ideal Self Affiliation Order Deference

Mean .606 .528 .602 .458 .441

SD .198 .209 .201 .247 .263

Min .010 .010 .010 .010 .050

(21)

Table 10. Estimated means and standard deviations of all estimated parameters for all five scales.

Scale

item Intraception Ideal Self Affiliation Order Deference α 2.02 (0.023) 2.10 (0.074) 2.07 (0.053) 2.04 (0.033) 2.27 (0.185) β 1 6.45 (1.215) 9.68 (3.273) 6.88 (0.968) 1.12 (0.378) 2.42 (0.332) 2 6.50 (1.001) 1.37 (0.257) 24.01 (6.884) 1.14 (0.263) 2.71 (0.625) 3 0.69 (0.332) 9.11 (5.300) 3.02 (3.731) 1.20 (0.265) 2.55 (0.857) 4 6.28 (1.197) 4.48 (5.629) 9.86 (4.483) 1.09 (0.277) 5.79 (2.155) 5 15.64 (5.344) 10.83 (3.375) 6.46 (4.009) 11.21 (4.101) 2.26 (1.062) 6 5.24 (0.643) 1.36 (1.475) 5.10 (0.824) 13.07 (10.717) 1.89 (0.418) 7 5.86 (0.653) 5.95 (2.426) 4.65 (0.908) 0.82 (0.330) 2.24 (4.581) 8 5.83 (0.938) 5.58 (0.820) 5.69 (0.876) 1.05 (0.357) 2.72 (0.447) 9 7.08 (0.653) 8.29 (5.592) 5.00 (0.854) 2.50 (0.508) 3.26 (1.233) 10 7.03 (1.732) 7.28 (1.931) 7.43 (2.160) 2.93 (0.582) 2.45 (0.322) γ 1 3.82 (0.589) 5.72 (1.776) 2.98 (0.348) 4.96 (4.305) 14.85 (2.675) 2 3.21 (0.415) 3.02 (1.466) 11.25 (3.000) 4.09 (2.465) 18.86 (6.010) 3 5.32 (6.380) 5.71 (2.726) 2.78 (3.710) 3.02 (1.894) 1.01 (0.148) 4 2.99 (0.461) 5.53 (5.578) 4.72 (1.885) 9.05 (4.333) 3.52 (1.151) 5 6.92 (2.216) 7.59 (2.186) 3.57 (1.910) 5.64 (1.879) 9.22 (7.881) 6 2.35 (0.238) 4.92 (5.578) 2.70 (0.332) 7.48 (4.655) 8.25 (3.728) 7 2.81 (0.275) 3.99 (1.508) 2.02 (0.251) 5.66 (6.532) 0.88 (0.306) 8 2.72 (0.372) 2.96 (0.357) 2.64 (0.318) 13.35 (10.005) 13.57 (3.015) 9 5.06 (0.997) 4.94 (3.239) 2.46 (0.311) 10.81 (3.303) 1.14 (0.160) 10 3.55 (0.734) 4.84 (1.072) 3.95 (0.967) 19.71 (5.235) 16.63 (3.084) ω 1 0.28 (0.090) 0.59 (0.089) 0.20 (0.058) 0.82 (0.123) 0.44 (0.044) 2 0.10 (0.057) 0.54 (0.192) 0.58 (0.058) 0.73 (0.111) 0.79 (0.039) 3 0.92 (0.072) 0.75 (0.134) 0.95 (0.068) 0.63 (0.174) 0.45 (0.194) 4 0.23 (0.075) 0.88 (0.135) 0.39 (0.107) 0.81 (0.057) 0.43 (0.134) 5 0.65 (0.064) 0.65 (0.085) 0.79 (0.123) 0.57 (0.082) 0.57 (0.300) 6 0.06 (0.043) 0.91 (0.106) 0.18 (0.082) 0.86 (0.078) 0.45 (0.132) 7 0.11 (0.050) 0.56 (0.176) 0.27 (0.093) 0.89 (0.088) 0.64 (0.251) 8 0.27 (0.067) 0.13 (0.064) 0.22 (0.069) 0.87 (0.067) 0.55 (0.050) 9 0.44 (0.096) 0.80 (0.111) 0.24 (0.083) 0.61 (0.063) 0.44 (0.171) 10 0.48 (0.088) 0.42 (0.102) 0.43 (0.104) 0.62 (0.044) 0.34 (0.043)

(22)

Table 11. Means and standard deviations of estimated parameters for the two items that are the

least biased and the two items that are the most biased.

Estimated parameters

Bias Scale Item Adjective ω α β γ

Small Affiliation 3 Balanced 0.95 (0.068) 2.07 (0.053) 3.02 (3.731) 2.78 (3.710) Intraception 3 Logical 0.92 (0.072) 2.02 (0.023) 0.69 (0.332) 5.32 (6.380) Large Intraception 6 Righteous 0.06 (0.043) 2.02 (0.023) 5.24 (0.643) 2.35 (0.238) Intraception 2 Intelligent 0.10 (0.057) 2.02 (0.023) 6.50 (1.001) 3.21 (0.415)

Figure 6. Distribution plots for the most unbiased items in green and for the most biased items in red.

As the Bayesian mixture model and standard one-factor model use different goodness of fit measures the two models will be compared using visual techniques. Histograms of the observed data are plotted with the expected curve of the model based on the estimated parameters. These histograms will be plotted side by side for the two models. A second comparison will be done using qq-plots. For the standard factor analysis the data are plotted against the expected normal distribution and for the mixture model the data are plotted against the mixture model based on the estimated parameters. Should the new model fit the data better than the standard maximum likelihood one-factor model the curve of the model should follow the distribution of the data more closely than the curve of the factor analysis will and the dots in the qq-plot should stay closer to the

(23)

diagonal line than they do for the standard factor analysis. All of these plots can be found in Appendix D.

It is clearly visible that while the mixture model does not always perfectly fit the data, its curve follows the distribution of the data much closer than the curve of the standard one-factor model does. The latter seems to only follow the distribution well, as expected, when it appears normal. The qq-plots show that the standard one-factor model fits the data rather well in the middle but under- and overestimates it quite often around the boundaries. The mixture model handles the data better, especially around the boundaries. These plots show that the mixture model that allows heteroscedastic data is a better suited option for line segments than a standard one-factor model using maximum likelihood estimation.

Both the simulation study and the application study show that the mixture model fits the data better than the standard one-factor model. First, the estimated parameters were close to their true values for most of the parameters, resulting in small RMSE values, indicating good estimation. Even the parameters with poorer recovery showed acceptable estimations. Second, the curve of the mixture model follows the distribution of the data much better than the curve of the standard one-factor model does. Therefore we introduce this model as the new method to fit a one-one-factor model on data measured using Visual Analogue Scales.

5. General Discussion

The goal of this study was to find a better method than the standard maximum likelihood factor model to fit a factor model to (heteroscedastic) line segment data. First, a data-based approach was tested. The arcsine transformation was applied to data on a [0,1] interval. It did not provide a better fit than the original data when a maximum likelihood one-factor model was fitted. Therefore the data-based approach was rejected and the next attempts to finding a better method focused on the creation of new heteroscedastic one-factor models.

Second, the heteroscedastic model presented in Hessen and Dolan (2007) and Molenaar et al (2010) was extended by changing the residual variance from a log linear relation between the error and the common factor to a beta function. While the beta function seemed promising for mild departures from normality the model did not work well on all heteroscedastic data.

Therefore last a mixture model was created and tested. The mixture model was a completely new and different model compared to the previous heteroscedastic models. It has one common factor that follows the beta function with equal shape parameters, creating a symmetrical

(24)

distribution, and an item-specific component. This component ensures that the model can handle heterogeneity in the data. This model is preferred over the standard one-factor model as it fits the data better. Therefore the conclusion of this study is that the beta mixture model should be preferred on line segment data over a maximum likelihood one-factor model.

However, there are some limitations to this study that should be noted. First of all this study only focuses on one-factor models. While it is an important model due to its application as an item response model for continuous data (Mellenbergh, 1994), many psychological variables are multidimensional. Therefore this model can only be used in very specific cases. However, it should not be too intricate to adjust this beta mixture model in such a way that it can fit multiple factors.

Second, because the mixture model is a Bayesian model and not that many goodness of fit measures have been developed for these models the output is limited compared to the standard one-factor model. For instance, the only goodness of fit measure provided by OpenBUGS is the deviance information criterion (DIC), while most factor analysis commands in R provide at least the Chi Square and RMSEA. Not only did it make it difficult to directly compare the new model to the old, it is also limiting when describing the fit of the model when applying it to data. However, in the former case the histograms and qq-plots provided in Appendix D are convincing enough that a lack of these goodness of fit measures is not a problem. On the other hand, the latter case still remains an issue. When more goodness of fit measures have been developed for Bayesian models the mixture model will be compared to the standard one-factor model more adequately.

Third, this lack of objective measures of goodness of fit results in subjectivity in the interpretation of the results. While the authors deemed the histograms and qq-plots convincing enough that the mixture model is a better suited method to fit a one-factor model on heteroscedastic data measured using line segments than the maximum likelihood one-factor model readers might disagree. More objective results might still be debatable, but give more certainty than the measures used in this study. However, in the meantime we are convinced enough by these results. Perhaps in the future the model can be tested more extensively.

Finally, the model was only tested on 15 items per scale in the simulation study and 10 items per scale in the application study. It is therefore the question how well the model will perform with different amounts of items. While there is no reason to believe that this should pose any problems in the future, it remains a possible limitation. Future applications of the model could indicate its capability when scales with different amounts of items are used.

In conclusion, this study presents a new heteroscedastic one-factor model that is especially suited for line segments on the interval [0,1]. We recommend this model instead of the standard one-factor model when fitting such a model to Visual Analogue Scales as it can handle heterogeneity in the data better than the latter, especially when the data is extreme.

(25)

6. Literature

Boker, S., Neale, M., Maes, H., Wilde, M., Spiegel M., Brick, T., Spies, J., Estabrook, R., Kenny, S., Bates, T., Mehta, P., & Fox, J. (2011). OpenMx: an open source extended structural equation modeling. Framework.Psychometrika, 76, 2, 306-317.

Bollen, K. A. (1996). A limited-information estimator for LISREL models with or without heteroscedastic errors. In G. A. Marcoulides & R. E. Schumacker (Eds.), Advanced structural

equation modeling: Issues and techniques (pp. 227-241). Mahwah, NJ: Erlbaum.

Digman, J. M. (1990). Personality structure: emergence of the five-factor model. Annual Review of

Psychology, 41, 1, 417-440.

Haviland, M. G., Warren, W. L., & Riggs, M. L. (2000). An observer scale to measure alexithymia. Psychosomatics, 41, 5, 385–392.

Hawker, G. A., Mian, S., Kendzerska, T., & French, M. (2011). Measures of adult pain: Visual Analogue Scale for pain (VAS Pain), Numeric Rating Scale for pain (NRS Pain), Mcgill Pain Questionnaire (MPQ), Short-Form McGill Pain Questionnaire (SF-MPQ), Chronic Pain Grade Scale (CPGS), Short Form-36 Bodily Pain Scale (SF-36 BPS), and measure of Intermittent and Constant Osteoarthritis Pain (ICOAP). Arthritis Care and Research, 63, 11, 240-252.

Hessen, D. J., & Dolan, C. V. (2007). Heteroscedastic one-factor models and marginal maximum likelihood estimation. The British Journal of Mathematical and Statistical Psychology, 62, 1, 57-77.

Hjermstad, M. J., Fayers, P. M., Haugen, D. F., Caraceni, A., Hanks, G. W., Loge, J. H., Fainsinger, R., Aass, N., & Kaasa, S. (2011). Studies comparing Numerical Rating Scales, Verbal Rating Scales, and Visual Analogue Scales for assessment of pain intensity in adults: a systematic literature review. Journal of Pain and Symptom Management, 41, 6, 1073-1093.

Goldberg, L. R. (1990). An alternative “description of personality”: the big-five factor structure.

Journal of Personality and Social Psychology, 59, 6, 1216.

Lawley, D. N. (1943). The application of the maximum likelihood method to factor analysis. British

Journal of Psychology, 33, 172-175.

Lawley, D. N., & Maxwell, A. E. (1971).Factor analysis as a statistical method. London: Butterworth. Lewin-Koh, S., & Amemiya, Y. (2003). Heteroscedastic factor analysis. Biometrika, 90, 85-97.

Meijer, E., & Mooijaart, A. (1996). Factor analysis with heteroscedastic errors. British Journal of

Mathematical and Statistical Psychology, 49, 189-202.

Mellenbergh, G. J. (1994). A unidimensional latent trait model for continuous item responses.

Multivariate Behavioral Research, 29, 3, 223-236.

Molenaar, D., Dolan, C. V., & Verhelst, N. D. (2010). Testing and modelling non-normality within the one-factor model. British Journal of Mathematical and Statistical Psychology, 63, 293-317.

(26)

Nelson, C. B., Rehm, J., Bedirhan, T., Grant, B., & Chatterji, S. (1999). Factor structures for DSM-IV substance disorder criteria endorsed by alcohol, cannabis, cocaine and opiate users: results from the WHO reliability and validity study. Addiction, 94, 6, 843-855.

Novick, M. R., & Jackson, P. H. (1974). Statistical methods for educational and psychological research. New York: McGraw-Hill.

Price, D. D., McGrath, P. A., Rafii, A., & Buckingham, B (1983). The validation of Visual Analogue Scales as ratio scale measures for chronic and experimental pain. Pain, 17, 1, 45-56.

Samejima, F. (1973). Homogeneous case of the continuous response model. Psychometrika, 38, 2, 203-219.

Taylor, G. J. (1984). Alexithymia: concept, measurement, and implications for treatment. The

American Journal of Psychiatry, 141, 6, 725-732.

Thomas, A., O’Hara, B., Ligges, U., & Sturtz, S. (2006). Making BUGS open. R News, 6, 12-17.

Vorst, H. C. M., & Bermond, B. (2011). Validity and reliability of the Bermond-Vorst Alexithymia Questionnaire. Personality and Individual Differences, 30, 3, 413-434.

(27)

Appendix A. R-code of the beta variance model

## simulate data Npeople=500 #people nit=15 #items

upsi=matrix(rep(c(.04,.02,.06,.03,.05),3),1,nit) #intercepts

lam=matrix(rep(c(0.90,0.93,0.92,0.94,0.95),3),1,nit) #factor loadings

indieta=matrix(round(rnorm(Npeople,mean=0.5,sd=0.1),4)) #normally distributed eta indiverr=matrix(,Npeople,nit) a=rep(2,15) b=rep(7,15) b0=rep(0.09,nit) for(p in 1:Npeople){ for(i in 1:nit){ indiverr[p,i]=rnorm(1,rep(0,nit),sqrt(b0[i]*indieta[p]^(a[i]-1)*(1-indieta[p])^(b[i]-1)))}} #normally distributed error with variance beta distributed

upsilon=matrix(rep(upsi,Npeople),Npeople,nit,byrow=T) data=upsilon+indieta%*%lam+indiverr

## use openMX to fit the model

# Following are prefixed, if changed, R-code below should be changed accordingly nfac=1 #number of factors

nq=30 #number of quadrature points # Setting things up library(OpenMx) nms=paste("y",1:nit) dimnames(data)[[2]]<-nms library(statmod) gauss.quad(nq,kind="hermite")->Quad Quad$weights/sqrt(pi)->W Quad$nodes*sqrt(2)->N

# Building up the OpenMx (sub)model(s) GenClass <- mxModel("GenClass",

(28)

mxMatrix( type="Full", nrow=nit, ncol=nfac, values=0, free=T,name="lamb" ), mxMatrix( type="Full", nrow=1, ncol=nit, free=T, values=0, name="b0"),

mxMatrix( type="Full", nrow=1, ncol=nit, free=T, values=2, name="a",lbound=0.01 ), mxMatrix( type="Full", nrow=1, ncol=nit, free=T, values=2, name="b",lbound=0.01 ), mxMatrix( type="Full", nrow=1, ncol=nit, values=0, free=T, name="upsil" ),

mxMatrix( type="Full", nrow=1, ncol=1, values=0, free=F, name="zeta"), mxMatrix( type="Full", nrow=1, ncol=1, values=0.01, free=F, name="varEta" ), mxMatrix( type="Full", nrow=1, ncol=1, values=0.5, free=F, name="meanEta" ))

createNewModel <- function(index, prefix, model) { modelname <- paste(prefix, index, sep='')

matNtmp <- mxMatrix(name="Ntmp",nrow=1,ncol=1,free=F,value=N[index])

matN <- mxAlgebra(name="N",expression=sqrt(container.omega2)*Ntmp+container.kappa) model <- mxModel(model, matN,matNtmp)

model <- mxRename(model, modelname) return(model)}

template <- mxModel("stClass",

mxMatrix(type="Full",nrow=nfac,ncol=nit,values=1,free=F,name="one"),

mxAlgebra( expression= GenClass.upsil + t(GenClass.lamb%x%N), name="condMean" ), mxAlgebra( expression=

(vec2diag(exp(GenClass.b0))%*%vec2diag(N%^%(GenClass.a-one))%*%vec2diag((1-N)%^%(GenClass.b-one))), name="condCov" ), mxFIMLObjective(covariance='condCov', means='condMean',vector=T,dimnames=nms)) topModel <- mxModel("container", mxMatrix(type="Full",nrow=nq,ncol=1,values=W,free=F,name="W"), mxMatrix(type="Full",nrow=nq,ncol=1,values=N,free=F,name="N"), mxAlgebra(expression=cbind(class1.objective,class2.objective,class3.objective,class4.objective,class5 .objective,class6.objective,class7.objective,class8.objective,class9.objective,class10.objective, class11.objective,class12.objective,class13.objective,class14.objective,class15.objective,class 16.objective,class17.objective,class18.objective,class19.objective,class20.objective,class21.o bjective,class22.objective,class23.objective,class24.objective,class25.objective,class26.object ive,class27.objective,class28.objective,class29.objective,class30.objective,class31.objective,cl ass32.objective,class33.objective,class34.objective,class35.objective,class36.objective,class3

(29)

7.objective,class38.objective,class39.objective,class40.objective,class41.objective,class42.obj ective,class43.objective,class44.objective,class45.objective,class46.objective,class47.objectiv e,class48.objective,class49.objective,class50.objective),name="LL"), mxMatrix(type="Full",nrow=nq,ncol=1,values=1,free=F,name="one"), mxAlgebra(expression=-log(2/(GenClass.zeta+1)-1),name="zetaStar"), mxAlgebra(expression= 2*sum(log(LL%*%((2%x%W)*(one/(one+((one+one)/((GenClass.zeta%x%one+one))-one)^(1.71%x%N)))))),name="defLL"), mxData(observed=data, type="raw"),

mxAlgebra(expression= GenClass.varEta / (1-2* zetaStar^2/(1+zetaStar^2)/pi), name="omega2"),

mxAlgebra(expression=GenClass.meanEta-sqrt(omega2*2/pi)*zetaStar/sqrt(zetaStar^2+1), name="kappa"),

mxAlgebraObjective("defLL"))

submodels <- lapply(1:nq, createNewModel, 'class', template) names(submodels) <- imxExtractNames(submodels)

topModel@submodels <- submodels topModel<-mxModel(topModel, GenClass)

# Fit the model and summarize results modelResults <- mxRun(topModel)

#save results of modelfit: u=summary(modelResults) u

(30)

Appendix B. R-code mixture model ## simulate data N=700 nit=15 x=3 eta=matrix(rbeta(N,x,x)) lamb=rep(.2,nit) b1=rep(1,nit) b2=rep(10,nit) err=y=matrix(,N,nit)

for(i in 1:nit) err[,i]=rbeta(N,b1,b2)

# the mixture for(p in 1:N) for(i in 1:nit){{ prob=sample(c(0,1),1,prob=c(1-lamb[i],lamb[i])) if(prob==0) y[p,i]=round(err[p,i],2) if(prob==1) y[p,i]=round(eta[p],2) if(y[p,i]==0.00) y[p,i]=0.01 if(y[p,i]==1.00) y[p,i]=0.99}} dat=list(N=N,nit=nit,y=y) ## starting values b1.init=matrix(c(2,NA),nit,2,T) b2.init=matrix(c(2,NA),nit,2,T) lamb.init=rep(0.5,nit) mode.init=matrix(sample(c(0,1),N*nit,T),N,nit) init=list(list(lamb=lamb.init,b1=b1.init,b2=b2.init,mode=mode.init,b.eta=2)) res=bugs(dat, init, c("b1","b2","lamb","b.eta"), 10000, model.file="directory", n.chains=1, n.burnin=5000,debug=F,codaPkg=F)

(31)

Appendix C: OpenBUGS-code mixture model model{ #data for(p in 1:N){ for(i in 1:nit){ mode[p,i]~dbern(lamb[i]) mode2[p,i]<-mode[p,i]+1 y[p,i]~dbeta(b1[i,mode2[p,i]],b2[i,mode2[p,i]]) eta[p,i]~dbeta(b1[i,2],b2[i,2]) }} #parameters for(i in 1:nit){ lamb[i]~dbeta(1,1) b1[i,1]~dgamma(.1,.1)I(.01,) b1[i,2]<-b.eta b2[i,2] <- b.eta b2[i,1]~dgamma(.1,.1)I(.01,) } #prior eta b.eta~dgamma(.1,.1)I(2,) }

(32)

Appendix D: parameter recovery of beta variance model

Table 12.

True values, mean observed values and RMSE values for all other parameters

υ λ

N items True Obs (SD) RMSE True Obs (SD) RMSE

300 1, 6, 11 .04 .12 (.001) .116 .90 .79 (.002) .148 2, 7, 12 .02 .10 (.002) .120 .93 .82 (.003) .151 3, 8, 13 .06 .14 (.001) .119 .92 .81 (.002) .151 4, 9, 14 .03 .115 (.001) .122 .94 .83 (.002) .155 5, 10, 15 .05 .14 (.001) .123 .95 .84 (.001) .155 500 1, 6, 11 .04 .09 (.001) .075 .90 .83 (.002) .129 2, 7, 12 .02 .07 (.001) .078 .93 .86 (.001) .134 3, 8, 13 .06 .11 (.001) .077 .92 .85 (.002) .135 4, 9, 14 .03 .08 (.001) .079 .94 .86 (.002) .135 5, 10, 15 .05 .10 (.001) .079 .95 .87 (.001) .136 700 1, 6, 11 .04 .09 (.001) .072 .90 .82 (.001) .129 2, 7, 12 .02 .07 (.001) .075 .93 .85 (.001) .133 3, 8, 13 .06 .11 (.001) .074 .92 .84 (.001) .132 4, 9, 14 .03 .08 (.001) .076 .94 .86 (.001) .136 5, 10, 15 .05 .10 (.001) .077 .95 .87 (.002) .137

(33)

Appendix E: mixture model plots

Standard FA, item 1

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0

Mixture Model, item 1

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y In tr a c e p ti o n

(34)

Standard FA, item 2 y In tr a c e p ti o n 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5

y In tr a c e p ti o n 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5 -3 -2 -1 0 1 2 3 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(35)

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(36)

Standard FA, item 4 y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(37)

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(38)

Standard FA, item 6 y In tr a c e p ti o n 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5

y In tr a c e p ti o n 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5 -3 -2 -1 0 1 2 3 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(39)

Standard FA, item 7 y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5 -3 -2 -1 0 1 2 3 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(40)

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(41)

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 2 .5 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(42)

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y In tr a c e p ti o n 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

(43)

Standard FA, item 1 y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(44)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(45)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(46)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(47)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(48)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(49)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(50)

Standard FA, item 8 y Id e a l S e lf 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0

y Id e a l S e lf 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(51)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf

(52)

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 1.0 0 .0 0 .5 1 .0 1 .5 2 .0 -3 -2 -1 0 1 2 3 0 .0 0 .2 0 .4 0 .6 0 .8

y Id e a l S e lf 0.0 0.2 0.4 0.6 0.8 0 .0 0 .2 0 .4 0 .6 0 .8 1 .0

y Id e a l S e lf