g03ca:: Multivariate Methods (NAG Toolbox)

Let

p

variables,

x_{1}, x_{2}, \dots, x_{p}

, with variance-covariance matrix

Σ

be observed. The aim of factor analysis is to account for the covariances in these

p

variables in terms of a smaller number,

k

, of hypothetical variables, or factors,

f_{1}, f_{2}, \dots, f_{k}

. These are assumed to be independent and to have unit variance. The relationship between the observed variables and the factors is given by the model:

x_{i} = \sum_{j = 1}^{k} λ_{i j} f_{j} + e_{i}, i = 1, 2, \dots, p

where

λ_{i j}

, for

i = 1, 2, \dots, p

and

j = 1, 2, \dots, k

, are the factor loadings and

e_{i}

, for

i = 1, 2, \dots, p

, are independent random variables with variances

ψ_{i}

, for

i = 1, 2, \dots, p

. The

ψ_{i}

represent the unique component of the variation of each observed variable. The proportion of variation for each variable accounted for by the factors is known as the communality. For this function it is assumed that both the

k

factors and the

e_{i}

's follow independent Normal distributions.

The estimation of the arguments of the model,

Λ

and

Ψ

, by maximum likelihood is described by Lawley and Maxwell (1971). The log-likelihood is:

- \frac{1}{2} (n - 1) \log (|Σ|) - \frac{1}{2} (n - 1) trace (S, Σ^{- 1}) + constant,

where

n

is the number of observations,

S

is the sample variance-covariance matrix or, if weights are used,

S

is the weighted sample variance-covariance matrix and

n

is the effective number of observations, that is, the sum of the weights. The constant is independent of the arguments of the model. A two stage maximization is employed. It makes use of the function

F (Ψ)

, which is, up to a constant,

- 2 / (n - 1)

times the log-likelihood maximized over

Λ

. This is then minimized with respect to

Ψ

to give the estimates,

\hat{Ψ}

, of

Ψ

. The function

F (Ψ)

can be written as:

F (Ψ) = \sum_{j = k + 1}^{p} (θ_{j} - \log θ_{j}) - (p - k)

where values

θ_{j}

, for

j = 1, 2, \dots, p

are the eigenvalues of the matrix:

S^{*} = Ψ^{- 1 / 2} S Ψ^{- 1 / 2} .

The estimates

\hat{Λ}

, of

Λ

, are then given by scaling the eigenvectors of

S^{*}

, which are denoted by

V

\hat{Λ} = Ψ^{1 / 2} V {(Θ - I)}^{1 / 2} .

where

Θ

is the diagonal matrix with elements

θ_{i}

, and

I

is the identity matrix.

The minimization of

F (Ψ)

is performed using nag_opt_bounds_mod_deriv2_comp (e04lb) which uses a modified Newton algorithm. The computation of the Hessian matrix is described by Clark (1970). However, instead of using the eigenvalue decomposition of the matrix

S^{*}

as described above, the singular value decomposition of the matrix

R Ψ^{- 1 / 2}

is used, where

R

is obtained either from the

Q R

decomposition of the (scaled) mean centred data matrix or from the Cholesky decomposition of the correlation/covariance matrix. The function nag_opt_bounds_mod_deriv2_comp (e04lb) ensures that the values of

ψ_{i}

are greater than a given small positive quantity,

δ

, so that the communality is always less than one. This avoids the so called Heywood cases.

In addition to the values of

Λ

Ψ

and the communalities, nag_mv_factor (g03ca) returns the residual correlations, i.e., the off-diagonal elements of

C - (Λ Λ^{T} + Ψ)

where

C

is the sample correlation matrix. nag_mv_factor (g03ca) also returns the test statistic:

χ^{2} = [n - 1 - (2 p + 5) / 6 - 2 k / 3] F (\hat{Ψ})

which can be used to test the goodness-of-fit of the model (1), see Lawley and Maxwell (1971) and Morrison (1967).

References

Parameters

Compulsory Input Parameters

Optional Input Parameters

Output Parameters

Error Indicators and Warnings

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

Accuracy

Further Comments

The factor loadings may be orthogonally rotated by using nag_mv_rot_orthomax (g03ba) and factor score coefficients can be computed using nag_mv_factor_score (g03cc). The maximum likelihood estimators are invariant to a change in scale. This means that the results obtained will be the same (up to a scaling factor) if either the correlation matrix or the variance-covariance matrix is used. As the correlation matrix ensures that all values of

ψ_{i}

are between

0

and

1

it will lead to a more efficient optimization. In the situation when the data matrix is input the results are always computed for the correlation matrix and then scaled if the results for the covariance matrix are required. When you input the covariance/correlation matrix the input matrix itself is used and you are advised to input the correlation matrix rather than the covariance matrix.

Example

function g03ca_example


fprintf('g03ca example results\n\n');

matrix = 'C';
n = int64(211);
x = [1,     0.523, 0.395, 0.471, 0.346, 0.426, 0.576, 0.434, 0.639;
     0.523, 1,     0.479, 0.506, 0.418, 0.462, 0.547, 0.283, 0.645;
     0.395, 0.479, 1,     0.355, 0.270, 0.254, 0.452, 0.219, 0.504;
     0.471, 0.506, 0.355, 1,     0.691, 0.791, 0.443, 0.285, 0.505;
     0.346, 0.418, 0.270, 0.691, 1,     0.679, 0.383, 0.149, 0.409;
     0.426, 0.462, 0.254, 0.791, 0.679, 1,     0.372, 0.314, 0.472;
     0.576, 0.547, 0.452, 0.443, 0.383, 0.372, 1,     0.385, 0.68;
     0.434, 0.283, 0.219, 0.285, 0.149, 0.314, 0.385, 1,     0.47;
     0.639, 0.645, 0.504, 0.505, 0.409, 0.472, 0.680, 0.470, 1];
nvar = int64(size(x,1));
isx  = ones(nvar,1,'int64');
nfac = int64(3);
iop = [int64(1); -1;  500;  2;  5];

[e, stat, com, psi, res, fl, ifail] = ...
  g03ca( ...
	 matrix, n, x, nvar, isx, nfac, iop);

disp(' Eigenvalues');
fprintf('%12.4e%12.4e%12.4e\n',e);
fprintf('\n     Test Statistic = %6.3f\n', stat(2));
fprintf('                 df = %6.3f\n', stat(3));
fprintf(' Significance level = %6.3f\n', stat(4));
fprintf('\n Residuals\n\n');
l = 1;
for i = 1:nvar-1
  fprintf('%8.3f', res(l:(l+i-1)));
  fprintf('\n');
  l = l + i;
end
fprintf('\n\n Loadings, Communalities and PSI\n\n');
for i = 1:nvar
  fprintf('%8.3f', fl(i,1:nfac), com(i), psi(i));
  fprintf('\n');
end

g03ca example results

 Eigenvalues
  1.5968e+01  4.3577e+00  1.8474e+00
  1.1560e+00  1.1190e+00  1.0271e+00
  9.2574e-01  8.9508e-01  8.7710e-01

     Test Statistic =  7.149
                 df = 12.000
 Significance level =  0.848

 Residuals

   0.000
  -0.013   0.022
   0.011  -0.005   0.023
  -0.010  -0.019  -0.016   0.003
  -0.005   0.011  -0.012  -0.001  -0.001
   0.015  -0.022  -0.011   0.002   0.029  -0.012
  -0.001  -0.011   0.013   0.005  -0.006  -0.001   0.003
  -0.006   0.010  -0.005  -0.011   0.002   0.007   0.003  -0.001


 Loadings, Communalities and PSI

   0.664  -0.321   0.074   0.550   0.450
   0.689  -0.247  -0.193   0.573   0.427
   0.493  -0.302  -0.222   0.383   0.617
   0.837   0.292  -0.035   0.788   0.212
   0.705   0.315  -0.153   0.619   0.381
   0.819   0.377   0.105   0.823   0.177
   0.661  -0.396  -0.078   0.600   0.400
   0.458  -0.296   0.491   0.538   0.462
   0.766  -0.427  -0.012   0.769   0.231

$iprint$	controls iteration monitoring;
	if $iprint \leq 0$ , then there is no printing of information else if $iprint > 0$ , then information is printed at every iprint iterations. The information printed consists of the value of $F (Ψ)$ at that iteration, the number of evaluations of $F (Ψ)$ , the current estimates of the communalities and an indication of whether or not they are at the boundary.
$maxfun$	the maximum number of function evaluations.
$acc$	the required accuracy for the estimates of $ψ_{i}$ .
$eps$	a lower bound for the values of $ψ$ , see Description.

On entry,	$ldfl < nvar$ ,
or	$nvar < 2$ ,
or	$n \leq nvar$ ,
or	$nfac < 1$ ,
or	$nvar < nfac$ ,
or	$m < nvar$ ,
or	$matrix ='D'$ or $'S'$ and $ldx < n$ ,
or	$matrix ='C'$ and $ldx < m$ ,
or	$matrix \neq'D'$ , $'S'$ or $'C'$ ,
or	$matrix ='D'$ or $'S'$ and $weight \neq'U'$ or $'W'$ ,
or	$iop (1) \neq 0$ and $iop (3)$ is such that $maxfun < 1$ ,
or	$iop (1) \neq 0$ and $iop (4)$ is such that $acc \geq 1.0$ ,
or	$iop (1) \neq 0$ and $iop (4)$ is such that $acc < machine precision$ ,
or	$iop (1) \neq 0$ and $iop (5)$ is such that $eps \geq 1.0$ ,
or	$iop (1) \neq 0$ and $iop (5)$ is such that $eps < machine precision$ ,
or	$matrix ='C'$ and $lwk < (5 \times nvar \times nvar + 33 \times nvar - 4) / 2$ ,
or	$matrix ='D'$ or $'S'$ and $lwk < \max ((5 \times nvar \times nvar + 33 \times nvar - 4) / 2, n \times nvar + 7 \times nvar + nvar \times (nvar - 1) / 2)$ .

NAG Toolbox: nag_mv_factor (g03ca)

▸▿ Contents

Purpose

Syntax

Description

References

Parameters

Compulsory Input Parameters

Optional Input Parameters

Output Parameters

Error Indicators and Warnings

Accuracy

Further Comments

Example