diff --git a/.pre-commit-config.yaml b/.pre-commit-config.yaml index 46c507a..ce90584 100644 --- a/.pre-commit-config.yaml +++ b/.pre-commit-config.yaml @@ -3,7 +3,7 @@ repos: - repo: https://github.com/pre-commit/pre-commit-hooks - rev: v4.1.0 + rev: v4.4.0 hooks: - id: trailing-whitespace - id: end-of-file-fixer @@ -14,27 +14,27 @@ repos: - id: check-json - id: debug-statements - id: check-merge-conflict -- repo: https://gitlab.com/pycqa/flake8 - rev: 4.0.1 +- repo: https://github.com/pycqa/flake8 + rev: 6.0.0 hooks: - id: flake8 args: [--max-line-length=101] exclude: ^(factor_analyzer/__init__.py) - repo: https://github.com/PyCQA/isort - rev: 5.10.1 + rev: 5.12.0 hooks: - id: isort args: ["--profile", "black", "--filter-files"] - repo: https://github.com/ikamensh/flynt/ - rev: '0.76' + rev: '0.78' hooks: - id: flynt - repo: https://github.com/psf/black - rev: 22.3.0 + rev: 23.3.0 hooks: - id: black - repo: https://github.com/pycqa/pydocstyle - rev: 6.1.1 + rev: 6.3.0 hooks: - id: pydocstyle args: diff --git a/factor_analyzer/factor_analyzer.py b/factor_analyzer/factor_analyzer.py index 47110b3..2d344d0 100644 --- a/factor_analyzer/factor_analyzer.py +++ b/factor_analyzer/factor_analyzer.py @@ -8,6 +8,7 @@ """ import warnings +from typing import Tuple import numpy as np import pandas as pd @@ -969,3 +970,36 @@ def get_factor_variance(self): check_is_fitted(self, "loadings_") loadings = self.loadings_.copy() return self._get_factor_variance(loadings) + + def sufficiency(self, num_observations: int) -> Tuple[float, int, float]: + """ + Perform the sufficiency test. + + The test calculates statistics under the null hypothesis that + the selected number of factors is sufficient. + + Parameters + ---------- + num_observations: int + The number of observations in the input data that this factor + analyzer was fit using. + + Returns + ------- + statistic: float + The test statistic + degrees: int + The degrees of freedom + pvalue: float + The p-value of the test + """ + nvar = self.corr_.shape[0] + degrees = ((nvar - self.n_factors) ** 2 - nvar - self.n_factors) // 2 + obj = self._fit_ml_objective( + self.get_uniquenesses(), self.corr_, self.n_factors + ) + statistic = ( + num_observations - 1 - (2 * nvar + 5) / 6 - (2 * self.n_factors) / 3 + ) * obj + pvalue = chi2.sf(statistic, df=degrees) + return statistic, degrees, pvalue diff --git a/tests/expected/test01/sufficiency_ml_none_15_test01.csv b/tests/expected/test01/sufficiency_ml_none_15_test01.csv new file mode 100644 index 0000000..666d7b3 --- /dev/null +++ b/tests/expected/test01/sufficiency_ml_none_15_test01.csv @@ -0,0 +1,16 @@ +n_factors,statistic,df,pvalue +1,6158.69420371228,740,0 +2,3244.72903681528,701,4.44659081257122e-322 +3,1475.87556412021,663,8.80428369753958e-64 +4,1206.41930867432,626,3.36923786804043e-39 +5,960.29751580042,590,3.64376268923982e-20 +6,820.739852110472,555,1.36586692869634e-12 +7,720.71679931586,521,1.37467963009699e-08 +8,631.742866356578,488,1.14373979137026e-05 +9,558.423820827789,456,0.000715862260222117 +10,500.223805331245,425,0.00686044615605591 +11,446.535371500571,395,0.0373055704171928 +12,397.200698549131,366,0.125803730285426 +13,343.188514729509,338,0.411317090149109 +14,304.335701489558,311,0.595750888646432 +15,269.814900486032,285,0.732252710489368 diff --git a/tests/sufficiency_test.r b/tests/sufficiency_test.r new file mode 100644 index 0000000..75f44fc --- /dev/null +++ b/tests/sufficiency_test.r @@ -0,0 +1,15 @@ +data = read.csv("tests/data/test01.csv") + +factor_range = 1:15 +statistic = rep(0, length(factor_range)) +dof = rep(0, length(factor_range)) +pvalue = rep(0, length(factor_range)) +for (i in 1:length(factor_range)) { + res = factanal(data, factors=factor_range[i], rotation="none") + statistic[i] = res$STATISTIC + dof[i] = res$dof + pvalue[i] = res$PVAL +} + +f = data.frame(n_factors=factor_range, statistic=statistic, df=dof, pvalue=pvalue) +write.csv(f, "tests/expected/test01/sufficiency_ml_none_15_test01.csv", row.names=FALSE, quote=FALSE) \ No newline at end of file diff --git a/tests/test_factor_analyzer.py b/tests/test_factor_analyzer.py index 24cda70..726619d 100644 --- a/tests/test_factor_analyzer.py +++ b/tests/test_factor_analyzer.py @@ -24,7 +24,6 @@ def test_calculate_bartlett_sphericity(): # noqa: D103 - path = "tests/data/test01.csv" data = pd.read_csv(path) s, p = calculate_bartlett_sphericity(data.values) @@ -34,7 +33,6 @@ def test_calculate_bartlett_sphericity(): # noqa: D103 def test_calculate_kmo(): # noqa: D103 - path = "tests/data/test02.csv" data = pd.read_csv(path) @@ -83,7 +81,6 @@ def test_gridsearch(): # noqa: D103 class TestFactorAnalyzer: # noqa: D101 def test_analyze_weights(self): # noqa: D102 - data = pd.DataFrame( { "A": [2, 4, 5, 6, 8, 9], @@ -107,7 +104,6 @@ def test_analyze_weights(self): # noqa: D102 assert_array_almost_equal(expected_weights, fa.weights_) def test_analyze_impute_mean(self): # noqa: D102 - data = pd.DataFrame( { "A": [2, 4, 5, 6, 8, 9], @@ -126,7 +122,6 @@ def test_analyze_impute_mean(self): # noqa: D102 assert_array_almost_equal(fa.corr_, expected_corr) def test_analyze_impute_median(self): # noqa: D102 - data = pd.DataFrame( { "A": [2, 4, 5, 6, 8, 9], @@ -148,7 +143,6 @@ def test_analyze_impute_median(self): # noqa: D102 assert_array_almost_equal(expected_corr, fa.corr_) def test_analyze_impute_drop(self): # noqa: D102 - data = pd.DataFrame( { "A": [2, 4, 5, 6, 8, 9], @@ -183,7 +177,6 @@ def test_analyze_rotation_value_error(self): # noqa: D102 @raises(ValueError) def test_analyze_infinite(self): # noqa: D102 - data = pd.DataFrame( { "A": [1.0, 0.4, 0.5], @@ -215,7 +208,6 @@ def test_smc_is_r_squared(self): # noqa: D102 assert_array_almost_equal(smc_result, expected_r2) def test_factor_variance(self): # noqa: D102 - path = "tests/data/test01.csv" data = pd.read_csv(path) @@ -234,3 +226,25 @@ def test_factor_variance(self): # noqa: D102 proportional_variance = fa.get_factor_variance()[1] assert_array_almost_equal(proportional_variance_expected, proportional_variance) + + def test_sufficiency(self): + path = "tests/data/test01.csv" + data = pd.read_csv(path) + + # compute the sufficiency values for a number of factors + computed_values = [] + for n_factors in range(1, 16): + fa = FactorAnalyzer(n_factors=n_factors, rotation=None, method="ml") + fa.fit(data) + computed_values.append([n_factors, *fa.sufficiency(data.shape[0])]) + + # create a data frame from the computed values + df_computed = pd.DataFrame( + computed_values, columns=["n_factors", "statistic", "df", "pvalue"] + ) + + # check against the values we expect + df_expected = pd.read_csv( + "tests/expected/test01/sufficiency_ml_none_15_test01.csv" + ) + pd.testing.assert_frame_equal(df_computed, df_expected) diff --git a/tests/test_utils.py b/tests/test_utils.py index ccc0e26..2c1d982 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -10,8 +10,7 @@ import numpy as np import pandas as pd from nose.tools import eq_, raises -from numpy.testing import assert_array_equal -from pandas.util.testing import assert_almost_equal +from numpy.testing import assert_array_equal, assert_almost_equal from factor_analyzer.utils import ( commutation_matrix, @@ -193,9 +192,9 @@ def test_partial_correlations(): # noqa: D103 data = pd.DataFrame([[12, 14, 15], [24, 12, 52], [35, 12, 41], [23, 12, 42]]) expected = [ - [1.0, -0.730955, -0.50616], - [-0.730955, 1.0, -0.928701], - [-0.50616, -0.928701, 1.0], + [1.0, -0.7309547, -0.50616], + [-0.7309547, 1.0, -0.9287013], + [-0.50616, -0.9287013, 1.0], ] expected = pd.DataFrame(expected, columns=[0, 1, 2], index=[0, 1, 2])