Equamax rotation giving similar results to quartimax instead

[nachomaiz]

Describe the bug Using equamax rotation gives the same results as quartimax to around 12 decimal points.

Quartimax results match with SPSS but not the equamax ones.

To Reproduce

import numpy as np
from factor_analyzer import FactorAnalyzer

data = np.random.randn(10,10)

equamax = FactorAnalyzer(n_factors=6, rotation='equamax', method='principal')
equamax.fit(data)
quartimax = FactorAnalyzer(n_factors=6, rotation='quartimax', method='principal')
quartimax.fit(data)

DECIMALS = 12

equal_loadings = np.array_equal(equamax.loadings_.round(DECIMALS),quartimax.loadings_.round(DECIMALS))
equal_rot_matrix = np.array_equal(equamax.rotation_matrix_.round(DECIMALS),quartimax.rotation_matrix_.round(DECIMALS))

print(equal_loadings, equal_rot_matrix)

This prints True True, but the loadings and rotation matrices should be different.

Expected behavior The equamax and quartimax results should be different from each other AND closely match with SPSS.

If you try the example above but replace 'quartimax' with 'varimax' the results are different.

Desktop (please complete the following information): Windows 10 x64, Python 3.9.7, factor_analyzer 0.3.2

[nachomaiz]

I've been looking at the code and, at a glance, it seems like it's running the right methods. No apparent typos when assigning the right objective for the Rotator class.

I also couldn't find a test for equamax in the tests folder.

I'm not too well-versed in the actual math behind it, but I guess if, according to the STATA manual equamax is some combination of varimax and quartimax, it should be different from both.

Maybe a dot product done backwards?

Hope this helps!

[jbiggsets]

Thanks! I'll look into this.

[nachomaiz]

Hey @jbiggsets,

I think I found the issue why I couldn't get equamax to match with SPSS outputs and why it was similar to quartimax.

It seems SPSS works out the kappa variable for the equamax rotation, while in every python implementation I've seen (scikit-learn, factor_analyzer, statsmodels, even R) kappa defaults to 0.

kappa in equamax is I believe n_factors / ( 2 * n_variables)

So, for the example above, SPSS would calculate a kappa of 6 / ( 2 * 10 ) = 0.3

Adding rotation_kwargs = {"kappa": (6 / 2 * data.shape[1])} to the FactorAnalyzer parameters in the example above solved the equality with quartimax and the mismatch with SPSS on my end.

I can see why one would want to leave kappa to be defined by the user since different rotations need different kappas.

Maybe a good compromise would be to make it more explicit in the documentation, or if the goal of the package is to replicate SPSS it might be worth it to hardcode the kappa calculation into _equamax_obj.

This was very helpful to figure it out: https://stackoverflow.com/questions/24300853/principal-component-analysis-with-equamax-rotation-in-r