From 3e3d653a23a85938a50e4e1e231aad0da363f195 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:03:24 +0100 Subject: [PATCH 1/9] BUG: Raise a clear error when a panel model has no residual degrees of freedom A model with as many parameters as observations once the effects are removed raised "ZeroDivisionError: division by zero" when the covariance was scaled by n / (n - k), and when the F-statistics were computed. The panels are small, but they are what is left of a panel after the missing values, the singletons and the absorbed effects are dropped: PooledOLS with 2 observations, a constant and a regressor, PanelOLS with entity or time effects for 2 entities with 2 periods and 2 regressors, BetweenOLS with 2 entities, a constant and a regressor, FirstDifferenceOLS with 1 entity of 3 periods and 2 regressors, and FamaMacBeth with a cross-section that has as many entities as regressors. The covariance now raises a ValueError that says that the model has no residual degrees of freedom, and how many observations, regressors and degrees of freedom of the effects it has. With debiased=False the estimates and the covariance are defined, the model F-statistic is an InvalidTestStatistic, and the summary can be printed. The pooled F-statistic of PanelOLS, which tests that the effects are zero, is an InapplicableTestStatistic if the effects do not use any degrees of freedom (a single entity or a single period with a constant, where the effect is the constant) and an InvalidTestStatistic if there are no residual degrees of freedom. The models with residual degrees of freedom do not change. The tests estimate each of the models with each covariance estimator, with and without debiased, and compare the slopes of the model with a single entity with the slopes of the pooled regression with a constant. --- linearmodels/panel/covariance.py | 7 + linearmodels/panel/model.py | 36 +++-- .../tests/panel/test_no_residual_df.py | 151 ++++++++++++++++++ 3 files changed, 184 insertions(+), 10 deletions(-) create mode 100644 linearmodels/tests/panel/test_no_residual_df.py diff --git a/linearmodels/panel/covariance.py b/linearmodels/panel/covariance.py index 880451024d..701234971c 100644 --- a/linearmodels/panel/covariance.py +++ b/linearmodels/panel/covariance.py @@ -105,6 +105,13 @@ def __init__( self._nobs_eff = self._nobs - extra_df if debiased: self._nobs_eff -= self._nvar + if self._nobs_eff <= 0: + regressors = f"{self._nvar} regressors and " if debiased else "" + raise ValueError( + "The covariance of the model cannot be estimated since the model has " + f"no residual degrees of freedom. It has {self._nobs} observations " + f"and uses {regressors}{extra_df} degrees of freedom for the effects." + ) self._scale = self._nobs / self._nobs_eff self._name = "Unadjusted" diff --git a/linearmodels/panel/model.py b/linearmodels/panel/model.py index a63ca58938..f3a9473fea 100644 --- a/linearmodels/panel/model.py +++ b/linearmodels/panel/model.py @@ -607,6 +607,10 @@ def _f_statistic( weps_const = y num_df = x.shape[1] name = "Model F-statistic (homoskedastic)" + if df_resid <= 0: + return InvalidTestStatistic( + "The model has no residual degrees of freedom", name=name + ) if self.has_constant: if num_df == 1: return InvalidTestStatistic("Model contains only a constant", name=name) @@ -2205,17 +2209,29 @@ def fit( df_num -= 1 weps_pooled = wy - wx @ _lstsq(wx, wy, rcond=None)[0] resid_ss_pooled = float(np.squeeze(weps_pooled.T @ weps_pooled)) - num = (resid_ss_pooled - resid_ss) / df_num - - denom = resid_ss / df_denom - stat = num / denom - f_pooled = WaldTestStatistic( - stat, - "Effects are zero", - df_num, - df_denom=df_denom, - name="Pooled F-statistic", + name = "Pooled F-statistic" + f_pooled: ( + WaldTestStatistic | InapplicableTestStatistic | InvalidTestStatistic ) + if df_num <= 0: + # A single entity (or period) with a constant: the effect is the constant + f_pooled = InapplicableTestStatistic( + reason="The effects do not use any degrees of freedom", name=name + ) + elif df_denom <= 0: + f_pooled = InvalidTestStatistic( + "The model has no residual degrees of freedom", name=name + ) + else: + num = (resid_ss_pooled - resid_ss) / df_num + denom = resid_ss / df_denom + f_pooled = WaldTestStatistic( + num / denom, + "Effects are zero", + df_num, + df_denom=df_denom, + name=name, + ) res.update(f_pooled=f_pooled) effects = DataFrame( eps_effects - eps, diff --git a/linearmodels/tests/panel/test_no_residual_df.py b/linearmodels/tests/panel/test_no_residual_df.py new file mode 100644 index 0000000000..abe21026ba --- /dev/null +++ b/linearmodels/tests/panel/test_no_residual_df.py @@ -0,0 +1,151 @@ +""" +Models without residual degrees of freedom + +A model that has as many parameters as observations, after the effects are removed, +raised the ZeroDivisionError of the scale of the covariance, n / (n - k), or of an +F-statistic. The covariance of such a model is not defined, and the error has to say +so. The F-statistic of the model is not defined either, which is the result if the +covariance is not debiased. The F-statistic that tests that the effects are zero is not +defined if the effects do not use any degree of freedom: a single entity with a +constant, or a single period with time effects. +""" + +import numpy as np +import pandas as pd +import pytest + +from linearmodels.panel import ( + BetweenOLS, + FamaMacBeth, + FirstDifferenceOLS, + PanelOLS, + PooledOLS, +) +from linearmodels.shared.exceptions import InferenceUnavailableWarning +from linearmodels.shared.hypotheses import ( + InapplicableTestStatistic, + InvalidTestStatistic, +) + +NO_DF = "no residual degrees of freedom" + + +def make_panel(nentity, nperiod, nvar, seed=0): + rg = np.random.default_rng(seed) + index = pd.MultiIndex.from_product( + [np.arange(nentity), np.arange(nperiod)], names=["entity", "time"] + ) + nobs = nentity * nperiod + x = pd.DataFrame( + rg.standard_normal((nobs, nvar)), + index=index, + columns=[f"x{i}" for i in range(nvar)], + ) + y = pd.Series(rg.standard_normal(nobs), index=index, name="y") + return y, x + + +def saturated_pooled(): + # 2 observations, a constant and a regressor + y, x = make_panel(1, 2, 1) + x.insert(0, "const", 1.0) + return y, x + + +def saturated_between(): + # 2 entities, a constant and a regressor in the regression of the means + y, x = make_panel(2, 3, 1) + x.insert(0, "const", 1.0) + return y, x + + +# Each model has as many parameters as observations once the effects are removed +SATURATED = { + "pooled": lambda: PooledOLS(*saturated_pooled()), + # 2 entities with 2 periods, 2 regressors and 2 effects + "entity-effects": lambda: PanelOLS(*make_panel(2, 2, 2), entity_effects=True), + "time-effects": lambda: PanelOLS(*make_panel(2, 2, 2), time_effects=True), + "between": lambda: BetweenOLS(*saturated_between()), + # 1 entity with 3 periods has 2 differences, and there are 2 regressors + "first-difference": lambda: FirstDifferenceOLS(*make_panel(1, 3, 2)), +} +# BetweenOLS does not have the kernel covariances +COVARIANCES = [ + (model, cov_type) + for model in SATURATED + for cov_type in (["robust"] if model == "between" else ["robust", "kernel", "ac"]) +] + + +@pytest.mark.parametrize("model", list(SATURATED)) +def test_saturated_model_raises_value_error(model): + mod = SATURATED[model]() + with pytest.raises(ValueError, match=NO_DF) as exc_info: + mod.fit() + assert exc_info.type is ValueError + + +@pytest.mark.parametrize(("model", "cov_type"), COVARIANCES) +def test_saturated_model_with_other_covariances(model, cov_type): + with pytest.raises(ValueError, match=NO_DF): + SATURATED[model]().fit(cov_type=cov_type) + + +@pytest.mark.parametrize("model", list(SATURATED)) +def test_saturated_model_that_is_not_debiased(model): + # The estimates and the covariance do not need the degrees of freedom. The + # F-statistic does, and it is not defined + res = SATURATED[model]().fit(debiased=False) + assert np.all(np.isfinite(res.params)) + assert isinstance(res.f_statistic, InvalidTestStatistic) + assert np.isnan(res.f_statistic.stat) + assert NO_DF in str(res.f_statistic) + assert str(res.summary) + + +def test_famamacbeth_with_as_many_observations_as_parameters(): + # One cross-section that has as many entities as regressors + y, x = make_panel(2, 1, 2) + with pytest.warns(InferenceUnavailableWarning): + mod = FamaMacBeth(y, x) + with pytest.raises(ValueError, match=NO_DF): + mod.fit() + + +def single_entity(const): + y, x = make_panel(1, 8, 2) + if const: + x.insert(0, "const", 1.0) + return y, x + + +def single_period(const): + y, x = make_panel(8, 1, 2) + if const: + x.insert(0, "const", 1.0) + return y, x + + +@pytest.mark.parametrize("const", [False, True]) +def test_single_entity_with_entity_effects(const): + # The effect is the constant: the slopes are those of the regression with a constant + y, x = single_entity(const) + res = PanelOLS(y, x, entity_effects=True).fit() + slopes = ["x0", "x1"] + pooled_x = x[slopes].copy() + pooled_x.insert(0, "const", 1.0) + pooled = PooledOLS(y, pooled_x).fit() + np.testing.assert_allclose(res.params[slopes], pooled.params[slopes]) + # The effects do not use any degree of freedom: there is nothing to test + assert isinstance(res.f_pooled, InapplicableTestStatistic) + assert "effects" in str(res.f_pooled) + assert np.isnan(res.f_pooled.stat) + assert str(res.summary) + + +@pytest.mark.parametrize("const", [False, True]) +def test_single_period_with_time_effects(const): + y, x = single_period(const) + res = PanelOLS(y, x, time_effects=True).fit() + assert isinstance(res.f_pooled, InapplicableTestStatistic) + assert np.all(np.isfinite(res.params)) From d3e82309e2da703de3c85e15a274ab804159b1a0 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:04:29 +0100 Subject: [PATCH 2/9] BUG: Check the rank of the means in BetweenOLS, and make df_model an integer BetweenOLS regresses the means of the entities. Regressors that are not collinear in the data can be collinear in the means: there can be fewer entities than regressors, and the dummy variables of the periods of a balanced panel all have the mean 1/T for every entity, so that they are collinear with the constant. The rank check of the model looks at the data, and the fit returned a minimum norm solution that is not an estimate of the parameters. The error came later, as numpy's "LinAlgError: Singular matrix" when the standard errors were used, and it did not say what was wrong. fit now checks the rank of the regression of the means, as the models check the rank of the data, and raises a ValueError that says that the regressors are collinear in the means of the entities, when it happens, and that check_rank=False skips the check. The other covariance estimators and the weights are not affected, and neither are the models that have the rank. The df_model of the results was the tuple (k,) and not the integer k that the other models return. The tests use a panel with more regressors than entities and the wage panel with the dummy variables of the years, which are the two cases, and check that the same models without the dummy variables are estimated with the degrees of freedom N - k. --- linearmodels/panel/model.py | 14 ++++- .../tests/panel/test_between_means.py | 60 +++++++++++++++++++ 2 files changed, 72 insertions(+), 2 deletions(-) create mode 100644 linearmodels/tests/panel/test_between_means.py diff --git a/linearmodels/panel/model.py b/linearmodels/panel/model.py index f3a9473fea..5b91ad5a0a 100644 --- a/linearmodels/panel/model.py +++ b/linearmodels/panel/model.py @@ -2437,10 +2437,20 @@ def fit( wx = root_w * x wy = root_w * y - params = _lstsq(wx, wy, rcond=None)[0] + params, _, rank, _ = _lstsq(wx, wy, rcond=None) + if self._check_rank and rank < x.shape[1]: + raise ValueError( + "The regressors are collinear in the means of the entities, which " + f"are the {y.shape[0]} observations of the regression, so that the " + "parameters are not identified. This happens if there are fewer " + "entities than regressors, or if a regressor is a combination of " + "the others only in the means, e.g., the dummy variables of the " + "periods and the constant in a balanced panel. If you wish to " + "proceed, you can set check_rank=False." + ) df_resid = y.shape[0] - x.shape[1] - df_model = (x.shape[1],) + df_model = x.shape[1] nobs = y.shape[0] cov_config = self._setup_clusters(cov_config) extra_df = 0 diff --git a/linearmodels/tests/panel/test_between_means.py b/linearmodels/tests/panel/test_between_means.py new file mode 100644 index 0000000000..de58f5a23c --- /dev/null +++ b/linearmodels/tests/panel/test_between_means.py @@ -0,0 +1,60 @@ +""" +BetweenOLS where the means of the regressors are collinear, and df_model + +The regression of BetweenOLS uses the means of the entities. The regressors can be +collinear in the means when they are not collinear in the data: there are fewer +entities than regressors, or the dummy variables of the periods are the constant 1/T in +the means of a balanced panel and so they are collinear with the constant. The +covariance raised numpy's "Singular matrix" error when it was used (and the parameters +were not identified), and the error has to say what is wrong. +""" + +import numpy as np +import pandas as pd +import pytest + +from linearmodels.datasets import wage_panel +from linearmodels.panel import BetweenOLS + + +def make_panel(nentity, nperiod, nvar, seed=0): + rg = np.random.default_rng(seed) + index = pd.MultiIndex.from_product( + [np.arange(nentity), np.arange(nperiod)], names=["entity", "time"] + ) + nobs = nentity * nperiod + x = pd.DataFrame( + rg.standard_normal((nobs, nvar)), + index=index, + columns=[f"x{i}" for i in range(nvar)], + ) + y = pd.Series(rg.standard_normal(nobs), index=index, name="y") + return y, x + + +def test_between_with_more_regressors_than_entities(): + y, x = make_panel(2, 3, 3) + with pytest.raises(ValueError, match="collinear in the means of the entities"): + BetweenOLS(y, x).fit() + # The rank is not checked if the user does not want it + res = BetweenOLS(y, x, check_rank=False).fit(debiased=False) + assert res.params.shape == (3,) + + +def test_between_with_time_dummies_in_a_balanced_panel(): + data = wage_panel.load() + year = pd.Categorical(data.year) + data = data.set_index(["nr", "year"]) + data["year"] = year + data["const"] = 1.0 + with pytest.raises(ValueError, match="collinear in the means of the entities"): + BetweenOLS(data.lwage, data[["const", "expersq", "union", "year"]]).fit() + res = BetweenOLS(data.lwage, data[["const", "expersq", "union"]]).fit() + assert res.df_resid == 545 - 3 + + +def test_between_df_model_is_the_number_of_regressors(): + y, x = make_panel(10, 3, 2) + res = BetweenOLS(y, x).fit() + assert res.df_model == 2 + assert isinstance(res.df_model, int) From dd85967e404fa1b85ac65b38c4651144cf4203f6 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:06:40 +0100 Subject: [PATCH 3/9] BUG: Report absorbed regressors for entity and time effects with one entity or period PanelOLS with entity_effects=True and time_effects=True for a panel that has a single entity (or a single period) raised pandas' "ValueError: No objects to concatenate". The two effects are removed by demeaning the data in the dimension with more groups, and then in the other dimension using the dummy variables of the dimension with fewer groups. With one entity or one period there are no dummy variables to remove, and concatenating the empty list failed. The data that have been demeaned in the first dimension are now returned when there are no dummy variables. The effects absorb all of the variation of the regressors, as they do for entity effects alone and a single entity, and the model raises the AbsorbingEffectError that says so, for each of the algorithms (the default, LSDV through use_lsdv, and low_memory). The tests check the error for panels of 1 entity and of 1 period, and that demeaning both dimensions of such a panel gives zeros. --- linearmodels/panel/data.py | 3 ++ .../panel/test_both_effects_degenerate.py | 49 +++++++++++++++++++ 2 files changed, 52 insertions(+) create mode 100644 linearmodels/tests/panel/test_both_effects_degenerate.py diff --git a/linearmodels/panel/data.py b/linearmodels/panel/data.py index aecca87421..a39f629729 100644 --- a/linearmodels/panel/data.py +++ b/linearmodels/panel/data.py @@ -442,6 +442,9 @@ def _demean_both(self, weights: PanelData | None) -> PanelData: e = self.demean(group, weights=weights) assert isinstance(e, PanelData) d = self.dummies(dummy, drop_first=True) + if d.shape[1] == 0: + # A single entity or period: there are no dummy variables to remove + return e d.index = e.index d_pd = PanelData(d).demean(group, weights=weights) d_arr = d_pd.values2d diff --git a/linearmodels/tests/panel/test_both_effects_degenerate.py b/linearmodels/tests/panel/test_both_effects_degenerate.py new file mode 100644 index 0000000000..3785765494 --- /dev/null +++ b/linearmodels/tests/panel/test_both_effects_degenerate.py @@ -0,0 +1,49 @@ +""" +Entity and time effects in a panel with a single entity or a single period + +The effects remove all of the variation of the regressors. The error has to say that the +effects absorbed the variables. The default algorithm that removes both effects +raised "No objects to concatenate" from pandas, since it uses the dummy variables of the +shorter dimension, and there are none. +""" + +import numpy as np +import pandas as pd +import pytest + +from linearmodels.panel import PanelOLS +from linearmodels.panel.data import PanelData +from linearmodels.panel.utility import AbsorbingEffectError + + +def make_panel(nentity, nperiod, nvar, seed=0): + rg = np.random.default_rng(seed) + index = pd.MultiIndex.from_product( + [np.arange(nentity), np.arange(nperiod)], names=["entity", "time"] + ) + nobs = nentity * nperiod + x = pd.DataFrame( + rg.standard_normal((nobs, nvar)), + index=index, + columns=[f"x{i}" for i in range(nvar)], + ) + y = pd.Series(rg.standard_normal(nobs), index=index, name="y") + return y, x + + +@pytest.mark.parametrize("low_memory", [None, False, True]) +@pytest.mark.parametrize("shape", [(1, 6), (6, 1)]) +def test_entity_and_time_effects_absorb_a_single_entity_or_period(shape, low_memory): + y, x = make_panel(*shape, 2) + mod = PanelOLS(y, x, entity_effects=True, time_effects=True) + with pytest.raises(AbsorbingEffectError, match="absorbed"): + mod.fit(low_memory=low_memory) + + +@pytest.mark.parametrize("shape", [(1, 6), (6, 1)]) +def test_demean_both_of_a_single_entity_or_period(shape): + # There is nothing to remove besides the means of the group with the observations + _, x = make_panel(*shape, 2) + demeaned = PanelData(x).demean("both") + assert demeaned.dataframe.shape == x.shape + np.testing.assert_allclose(demeaned.values2d, 0.0, atol=1e-14) From 8b13c269b7fca0348eec007663d0b63fccb09d74 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:08:12 +0100 Subject: [PATCH 4/9] BUG: Accept a NumPy Generator in generate_panel_data generate_panel_data(rng=numpy.random.default_rng(0)) raised "AttributeError: 'Generator' object has no attribute 'randint'". The function draws the categories of the effects and of the clusters with randint, which only a RandomState has, and the documentation only said RandomState. Generators are what NumPy recommends, and what the rest of the package and its tests use, so the function now draws the integers with integers for a Generator and with randint for a RandomState. The other draws (standard_normal, chisquare and choice) are methods of both. Calls with a RandomState, or without rng, give the same data as before. The existing test of the function listed random states as a parameter but never passed them to the function. It does now, with a Generator added, and a new test checks that a seed gives the same data, the ranges of the categories, and the values that a RandomState gives (taken from the function before the change). --- linearmodels/panel/utility.py | 23 ++++++++--- linearmodels/tests/panel/test_utility.py | 50 +++++++++++++++++++++++- 2 files changed, 65 insertions(+), 8 deletions(-) diff --git a/linearmodels/panel/utility.py b/linearmodels/panel/utility.py index 6f2e47cbd4..10a5613621 100644 --- a/linearmodels/panel/utility.py +++ b/linearmodels/panel/utility.py @@ -642,6 +642,17 @@ class PanelModelData(NamedTuple): clusters: DataFrame +def _integers( + rng: numpy.random.RandomState | numpy.random.Generator, + high: int, + size: tuple[int, ...], +) -> linearmodels.typing.data.IntArray: + """Random integers in [0, high) from a RandomState or a Generator""" + if isinstance(rng, numpy.random.Generator): + return rng.integers(0, high, size) + return rng.randint(0, high, size) + + def generate_panel_data( nentity: int = 971, ntime: int = 7, @@ -650,7 +661,7 @@ def generate_panel_data( missing: float = 0, other_effects: int = 2, ncats: int | list[int] = 4, - rng: numpy.random.RandomState | None = None, + rng: numpy.random.RandomState | numpy.random.Generator | None = None, ) -> PanelModelData: """ Simulate panel data for testing @@ -673,9 +684,9 @@ def generate_panel_data( The number of categories to use in other_effects and variance clusters. If list-like, then it must have as many elements as other_effects. - rng : RandomState - A NumPy RandomState instance. If not provided, one is initialized - using a fixed seed. + rng : RandomState or Generator + A NumPy RandomState or Generator instance. If not provided, a + RandomState is initialized using a fixed seed. Returns ------- @@ -721,11 +732,11 @@ def generate_panel_data( _c = [] for i in range(other_effects): nc = ncats[i] - _c.append(rng.randint(0, nc, (1, t, n))) + _c.append(_integers(rng, nc, (1, t, n))) c = np.concatenate(_c, 0) vcats = [f"varcat.{i}" for i in range(2)] - vc2 = np.ones((2, t, 1)) @ rng.randint(0, n // 2, (2, 1, n)) + vc2 = np.ones((2, t, 1)) @ _integers(rng, n // 2, (2, 1, n)) vc1 = vc2[[0]] if const: diff --git a/linearmodels/tests/panel/test_utility.py b/linearmodels/tests/panel/test_utility.py index 8613a473e8..826d139918 100644 --- a/linearmodels/tests/panel/test_utility.py +++ b/linearmodels/tests/panel/test_utility.py @@ -224,7 +224,9 @@ class SubArray(np.ndarray): @pytest.mark.parametrize("const", [True, False]) @pytest.mark.parametrize("other_effects", [0, 1, 2]) @pytest.mark.parametrize("cat_list", [True, False]) -@pytest.mark.parametrize("random_state", [None, np.random.RandomState(0)]) +@pytest.mark.parametrize( + "random_state", [None, np.random.RandomState(0), np.random.default_rng(0)] +) def test_generate_panel_data(missing, const, other_effects, cat_list, random_state): if cat_list: ncats: list[int] | int = [13] * other_effects @@ -232,7 +234,11 @@ def test_generate_panel_data(missing, const, other_effects, cat_list, random_sta ncats = 21 dataset = generate_panel_data( - missing=missing, const=const, other_effects=other_effects, ncats=ncats + missing=missing, + const=const, + other_effects=other_effects, + ncats=ncats, + rng=random_state, ) assert isinstance(dataset, PanelModelData) @@ -244,6 +250,46 @@ def test_generate_panel_data(missing, const, other_effects, cat_list, random_sta assert dataset.other_effects.shape == (dataset.data.shape[0], other_effects) +@pytest.mark.parametrize("make_rng", [np.random.RandomState, np.random.default_rng]) +def test_generate_panel_data_rng(make_rng): + # A Generator has integers and not randint. The same seed gives the same data + options = {"nentity": 30, "ntime": 5, "other_effects": 2, "ncats": [3, 4]} + first = generate_panel_data(rng=make_rng(1234), **options) + second = generate_panel_data(rng=make_rng(1234), **options) + other = generate_panel_data(rng=make_rng(4321), **options) + for result in (first, second, other): + assert result.data.shape == (150, 6) + assert result.other_effects.shape == (150, 2) + pd.testing.assert_frame_equal(first.data, second.data) + pd.testing.assert_frame_equal(first.other_effects, second.other_effects) + pd.testing.assert_frame_equal(first.clusters, second.clusters) + assert not first.data.equals(other.data) + # The categories of the effects are integers in the ranges + effects = first.other_effects + assert set(effects["cat.0"].unique()) <= set(range(3)) + assert set(effects["cat.1"].unique()) <= set(range(4)) + + +def test_generate_panel_data_random_state_is_unchanged(): + # The draws of a RandomState do not depend on Generators being accepted + result = generate_panel_data( + nentity=6, + ntime=3, + nexog=2, + other_effects=2, + ncats=[3, 5], + rng=np.random.RandomState(7), + ) + effects = result.other_effects.to_numpy().T + assert_array_equal( + effects, + [ + [1, 0, 1, 2, 1, 0, 1, 1, 0, 2, 0, 2, 2, 1, 1, 0, 1, 1], + [2, 1, 1, 0, 3, 1, 2, 1, 3, 0, 0, 2, 2, 0, 0, 1, 1, 2], + ], + ) + + def test_not_absorbed_const(): x = np.random.standard_normal((200, 3)) x[:, 0] = 0 From 7bb078dbfab290f8e87832660a995bedd743fe53 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:09:23 +0100 Subject: [PATCH 5/9] BUG: Return one row of clusters for each observation from generate_panel_data The clusters of generate_panel_data were the concatenation, along the rows, of a frame with the first cluster variable and of a frame with both cluster variables. The result had twice as many rows as the data, with an index that is not unique, and half of the rows had a missing value of the second variable. Clustered covariances need one row of clusters for each observation, and the frame could not be used with the data that the function returns. The clusters are now the frame with both variables, varcat.0 and varcat.1, for the observations of the data. The first variable is the one that the removed frame had, so that the values are those that a RandomState gives for the second half of the old frame, and the draws of the other variables are unchanged. The tests check the shape, the index, the missing values, the columns, the range of the values, that the clusters of an entity are the same in all of its periods, the values for a seeded RandomState (the values of the old frame), and that the clusters work in a clustered covariance of PanelOLS. --- linearmodels/panel/utility.py | 10 ++----- linearmodels/tests/panel/test_utility.py | 36 ++++++++++++++++++++++++ 2 files changed, 38 insertions(+), 8 deletions(-) diff --git a/linearmodels/panel/utility.py b/linearmodels/panel/utility.py index 10a5613621..8e0470aa34 100644 --- a/linearmodels/panel/utility.py +++ b/linearmodels/panel/utility.py @@ -737,7 +737,6 @@ def generate_panel_data( vcats = [f"varcat.{i}" for i in range(2)] vc2 = np.ones((2, t, 1)) @ _integers(rng, n // 2, (2, 1, n)) - vc1 = vc2[[0]] if const: x[0] = 1.0 @@ -771,14 +770,9 @@ def generate_panel_data( c, items=cats, major_axis=time, minor_axis=entities, swap=True ) other_eff = c_df.reindex(index) - vc1_df = panel_to_frame( - vc1, items=vcats[:1], major_axis=time, minor_axis=entities, swap=True - ) - vc1_df = vc1_df.reindex(index) - vc2_df = panel_to_frame( + clusters = panel_to_frame( vc2, items=vcats, major_axis=time, minor_axis=entities, swap=True ) - vc2_df = vc2_df.reindex(index) - clusters = concat([vc1_df, vc2_df], sort=False) + clusters = clusters.reindex(index) data = concat([y_df, x_df], axis=1, sort=False) return PanelModelData(data, w_df, other_eff, clusters) diff --git a/linearmodels/tests/panel/test_utility.py b/linearmodels/tests/panel/test_utility.py index 826d139918..723fcfa127 100644 --- a/linearmodels/tests/panel/test_utility.py +++ b/linearmodels/tests/panel/test_utility.py @@ -5,6 +5,7 @@ import pandas as pd import pytest +from linearmodels.panel import PanelOLS from linearmodels.panel.utility import ( AbsorbingEffectError, PanelModelData, @@ -290,6 +291,41 @@ def test_generate_panel_data_random_state_is_unchanged(): ) +def test_generate_panel_data_clusters_random_state_is_unchanged(): + # The clusters are the draws that the second column of the clusters had + result = generate_panel_data( + nentity=6, + ntime=3, + nexog=2, + other_effects=2, + ncats=[3, 5], + rng=np.random.RandomState(7), + ) + first_periods = result.clusters.iloc[::3] + assert_array_equal(first_periods["varcat.0"], [2, 2, 0, 2, 0, 0]) + assert_array_equal(first_periods["varcat.1"], [1, 0, 1, 0, 2, 2]) + + +@pytest.mark.parametrize("missing", [0, 0.1]) +def test_generate_panel_data_clusters(missing): + # The clusters have one row for each observation, and two variables + result = generate_panel_data(nentity=20, ntime=5, missing=missing) + assert result.clusters.shape == (result.data.shape[0], 2) + assert result.clusters.index.equals(result.data.index) + assert not result.clusters.isnull().any().any() + assert result.clusters.to_numpy().max() < 10 + assert list(result.clusters.columns) == ["varcat.0", "varcat.1"] + # They are the clusters of the entities: the same in all periods of an entity + assert (result.clusters.groupby(level=0).nunique() == 1).all().all() + # And can be used in a covariance + data = result.data.dropna() + clusters = result.clusters.loc[data.index] + res = PanelOLS(data["y"], data[["x0", "x1"]], entity_effects=True).fit( + cov_type="clustered", clusters=clusters + ) + assert np.all(np.isfinite(res.params)) + + def test_not_absorbed_const(): x = np.random.standard_normal((200, 3)) x[:, 0] = 0 From a68e4717d3548d6f11170e3d654e9c36f5ba2cbf Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:10:55 +0100 Subject: [PATCH 6/9] DOC: Describe the FamaMacBeth regressions as cross-sectional in the warnings The two warnings of FamaMacBeth about the regressions that cannot be estimated called them "time-series regressions" (and "time-series observation available to estimate cross-sectional regressions"). The model estimates one cross-sectional regression for each period, and then uses the time series of the estimates. The warnings now say how many cross-sectional regressions cannot be estimated due to deficient rank, and how many can be estimated when there are fewer than there are parameters, which is when inference is not available. The tests check the texts of the warnings. --- linearmodels/panel/model.py | 9 ++++----- linearmodels/tests/panel/test_fama_macbeth.py | 13 ++++++++++--- 2 files changed, 14 insertions(+), 8 deletions(-) diff --git a/linearmodels/panel/model.py b/linearmodels/panel/model.py index 5b91ad5a0a..73dfb42696 100644 --- a/linearmodels/panel/model.py +++ b/linearmodels/panel/model.py @@ -3337,17 +3337,16 @@ def _mr(ex: linearmodels.typing.data.Float64Array) -> int: if valid_blocks.sum() < exog.shape[1]: warnings.warn( - "The number of time-series observation available to estimate " - "cross-sectional\nregressions, {}, is less than the number of " - "parameters in the model. Parameter\ninference is not " - "available.".format(valid_blocks.sum()), + "The number of cross-sectional regressions that can be estimated, " + "{}, is less than the number of parameters in the model. Parameter " + "inference is not available.".format(valid_blocks.sum()), InferenceUnavailableWarning, stacklevel=3, ) elif valid_blocks.sum() < valid_blocks.shape[0]: warnings.warn( - "{} of the time-series regressions cannot be estimated due to " + "{} of the cross-sectional regressions cannot be estimated due to " "deficient rank.".format(valid_blocks.shape[0] - valid_blocks.sum()), MissingValueWarning, stacklevel=3, diff --git a/linearmodels/tests/panel/test_fama_macbeth.py b/linearmodels/tests/panel/test_fama_macbeth.py index 6da1ab917c..16ba6342f7 100644 --- a/linearmodels/tests/panel/test_fama_macbeth.py +++ b/linearmodels/tests/panel/test_fama_macbeth.py @@ -129,9 +129,14 @@ def test_block_size_warnings(): ) y = pd.DataFrame(y, index=idx, columns=["y"]) x = pd.DataFrame(x, index=idx, columns=["x1", "x2", "x3"]) - with pytest.warns(MissingValueWarning): + # The periods are cross-sectional regressions, and not time-series regressions + with pytest.warns( + MissingValueWarning, match=r"1 of the cross-sectional regressions" + ): FamaMacBeth(y.iloc[:11], x.iloc[:11]) - with pytest.warns(InferenceUnavailableWarning): + with pytest.warns( + InferenceUnavailableWarning, match=r"cross-sectional regressions that can be" + ): FamaMacBeth(y.iloc[::4], x.iloc[::4]) @@ -338,7 +343,9 @@ def test_inference_n40_t10(): def test_inference_dropped_periods(bandwidth, debiased): # T is the number of periods that are estimated, 21 of the 25 periods y, x = dropped_periods_panel() - with pytest.warns(MissingValueWarning, match=r"4 of the time-series regressions"): + with pytest.warns( + MissingValueWarning, match=r"4 of the cross-sectional regressions" + ): mod = FamaMacBeth(y, x) cov_type = "kernel" if bandwidth else "unadjusted" res = mod.fit( From cf1f7ed4635ccf07d51e47c6ef0aed18a9dc1334 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:13:27 +0100 Subject: [PATCH 7/9] BUG: Use the number of regressors as the rank when check_rank=False With check_rank=False the panel models skip the computation of the rank of exog and use the number of columns, which is the rank if the regressors are not collinear. They used self.exog.shape[1], and the shape of the PanelData is (variables, periods, entities), so that the rank was the number of periods. The rank is used to decide if the regressors contain a constant that is not a column, because the regressors with a constant have the rank of the regressors alone. The wrong rank gave the wrong answer: * Two dummy variables that add to one, in a panel with 5 periods, were not a constant. The R-squared of PooledOLS was uncentered, 0.2942 instead of 0.1877, that of BetweenOLS 0.5137 instead of 0.0709, and the F-statistics tested the model as if it had no constant. * Three regressors without a constant, in a panel with 2 periods, were a constant. The R-squared of PooledOLS was -0.8760 instead of 0.0304 and FirstDifferenceOLS raised "Constants are not allowed in first difference regressions". The models whose regressors have a constant column, which is found before the rank is used, and the default check_rank=True do not change. The tests compare the models with check_rank=False with those with check_rank=True, and the R-squared of PooledOLS with the R-squared written out in NumPy, for both panels. --- linearmodels/panel/model.py | 5 +- .../tests/panel/test_check_rank_false.py | 124 ++++++++++++++++++ 2 files changed, 127 insertions(+), 2 deletions(-) create mode 100644 linearmodels/tests/panel/test_check_rank_false.py diff --git a/linearmodels/panel/model.py b/linearmodels/panel/model.py index 73dfb42696..9b635f8a7f 100644 --- a/linearmodels/panel/model.py +++ b/linearmodels/panel/model.py @@ -523,9 +523,10 @@ def _grid_to_rows(self, grid: PanelData, labelled: bool) -> PanelData: return PanelData(rows) def _check_exog_rank(self) -> int: - if not self._check_rank: - return self.exog.shape[1] x = cast("linearmodels.typing.data.Float64Array", self.exog.values2d) + if not self._check_rank: + # The shape of the PanelData is (variables, periods, entities) + return x.shape[1] _, _, rank_of_x, _ = _lstsq(x, np.ones(x.shape[0])) if rank_of_x < x.shape[1]: raise ValueError( diff --git a/linearmodels/tests/panel/test_check_rank_false.py b/linearmodels/tests/panel/test_check_rank_false.py new file mode 100644 index 0000000000..4d35c85886 --- /dev/null +++ b/linearmodels/tests/panel/test_check_rank_false.py @@ -0,0 +1,124 @@ +""" +The constant is detected whether or not the rank of the regressors is checked + +check_rank=False skips the computation of the rank of the regressors, and the models +then use the number of regressors, which is the rank if they are not collinear. They +used the number of periods of the panel instead. The detection of a constant that is +not a column, e.g., dummy variables that add to one, compares the rank of the +regressors with the rank of the regressors and a constant, so that it decided that +there was no constant in a panel with 5 periods and 2 dummy variables, and that there +was a constant in a panel with 2 periods and 3 regressors. The R-squared and the +F-statistic then treated the models as if they had, or did not have, a constant. +""" + +import numpy as np +from numpy.testing import assert_allclose +import pandas as pd +import pytest + +from linearmodels.panel import ( + BetweenOLS, + FamaMacBeth, + FirstDifferenceOLS, + PanelOLS, + PooledOLS, + RandomEffects, +) + +MODELS = ["pooled", "between", "random-effects", "fama-macbeth", "panel-ols"] + +# Two cross-sections cannot give the standard errors of 3 parameters +pytestmark = pytest.mark.filterwarnings( + "ignore:The number of cross-sectional regressions:" + "linearmodels.shared.exceptions.InferenceUnavailableWarning" +) + + +def build(model, y, x, check_rank): + if model == "panel-ols": + return PanelOLS(y, x, entity_effects=True, check_rank=check_rank) + models = { + "pooled": PooledOLS, + "between": BetweenOLS, + "random-effects": RandomEffects, + "fama-macbeth": FamaMacBeth, + "first-difference": FirstDifferenceOLS, + } + return models[model](y, x, check_rank=check_rank) + + +def make_index(nentity, nperiod): + return pd.MultiIndex.from_product( + [np.arange(nentity), np.arange(nperiod)], names=["entity", "time"] + ) + + +def implicit_constant_panel(): + """Two dummy variables that add to one, and 5 periods""" + rg = np.random.default_rng(0) + index = make_index(40, 5) + group = rg.integers(0, 2, len(index)) + x = pd.DataFrame( + {"a": (group == 0) * 1.0, "b": (group == 1) * 1.0}, + index=index, + ) + y = pd.Series(rg.standard_normal(len(index)) + group, index=index, name="y") + return y, x + + +def no_constant_panel(): + """Three regressors with a mean of zero, and 2 periods""" + rg = np.random.default_rng(1) + index = make_index(40, 2) + x = pd.DataFrame( + rg.standard_normal((len(index), 3)), index=index, columns=["a", "b", "c"] + ) + y = pd.Series(rg.standard_normal(len(index)) + 1, index=index, name="y") + return y, x + + +def r_squared(y, x, centered): + yv, xv = y.to_numpy(), x.to_numpy() + resid = yv - xv @ np.linalg.lstsq(xv, yv, rcond=None)[0] + total = yv - yv.mean() if centered else yv + return 1 - resid @ resid / (total @ total) + + +@pytest.mark.parametrize("model", MODELS) +def test_implicit_constant(model): + y, x = implicit_constant_panel() + checked = build(model, y, x, True) + skipped = build(model, y, x, False) + assert checked.has_constant + assert skipped.has_constant + res_checked, res_skipped = checked.fit(), skipped.fit() + assert_allclose(res_skipped.params, res_checked.params) + assert_allclose(res_skipped.rsquared, res_checked.rsquared) + assert_allclose(res_skipped.f_statistic.stat, res_checked.f_statistic.stat) + + +@pytest.mark.parametrize("model", MODELS + ["first-difference"]) +def test_no_constant_with_more_regressors_than_periods(model): + y, x = no_constant_panel() + checked = build(model, y, x, True) + skipped = build(model, y, x, False) + assert not checked.has_constant + assert not skipped.has_constant + res_checked, res_skipped = checked.fit(), skipped.fit() + assert_allclose(res_skipped.params, res_checked.params) + assert_allclose(res_skipped.rsquared, res_checked.rsquared) + assert_allclose(res_skipped.f_statistic.stat, res_checked.f_statistic.stat) + + +@pytest.mark.parametrize("check_rank", [True, False]) +def test_pooled_r_squared_of_an_implicit_constant_is_centered(check_rank): + y, x = implicit_constant_panel() + res = PooledOLS(y, x, check_rank=check_rank).fit() + assert_allclose(res.rsquared, r_squared(y, x, True)) + + +@pytest.mark.parametrize("check_rank", [True, False]) +def test_pooled_r_squared_without_a_constant_is_not_centered(check_rank): + y, x = no_constant_panel() + res = PooledOLS(y, x, check_rank=check_rank).fit() + assert_allclose(res.rsquared, r_squared(y, x, False)) From 07e775aa73ea1a057dd533999a7a5d9ad69dc6d2 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:26:28 +0100 Subject: [PATCH 8/9] TST: Expect the constant error of FirstDifferenceOLS when the rank is not checked test_formulas_rank_check fits each model with check_rank=False to regressors that are not full rank, and FirstDifferenceOLS passed only because the rank that the models used when the rank was not checked was the number of periods of the panel and not the number of regressors. The function that looks for a constant that is not a column reads regressors that are not full rank as containing an implicit constant, and first differences do not allow a constant. With the number of regressors as the rank, FirstDifferenceOLS raises "Constants are not allowed in first difference regressions." for such regressors, as it does when the rank is computed. The test now expects the error for that model. The other models do not change: the regressors with check_rank=False are estimated as before. --- linearmodels/tests/panel/test_formula.py | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/linearmodels/tests/panel/test_formula.py b/linearmodels/tests/panel/test_formula.py index 7897ae6ec7..6573ab5636 100644 --- a/linearmodels/tests/panel/test_formula.py +++ b/linearmodels/tests/panel/test_formula.py @@ -334,6 +334,16 @@ def _mr(z): # FMB has tighter checks that cannot be disabled return + if model == FirstDifferenceOLS: + # The check for a constant reads regressors that are not full rank as having an + # implicit constant, which first differences do not allow. The rank that it + # used when the rank was not checked was the number of periods, not the rank + with pytest.raises(ValueError, match=r"Constants are not allowed"): + model.from_formula(formula, joined, check_rank=False) + with pytest.raises(ValueError, match=r"Constants are not allowed"): + model(y, x, check_rank=False) + return + mod = model.from_formula(formula, joined, check_rank=False) res = mod.fit() assert isinstance(res.summary.as_text(), str) From c3b49b6980f6f4ac74b870116816f05a382b88e2 Mon Sep 17 00:00:00 2001 From: Kevin Sheppard Date: Mon, 12 Oct 2026 20:42:27 +0100 Subject: [PATCH 9/9] BUG: Reject infinite weights in the panel models The panel models check that the weights are strictly positive, once the observations with missing weights have been dropped, and then scale them to have a mean of one. An infinite weight passed the check, and the scaling made it NaN and every other weight zero. The models then failed with errors from NumPy, SciPy and pandas that do not say what is wrong: "array must not contain infs or NaNs", "SVD did not converge", an AbsorbingEffectError for entity effects, and "operands could not be broadcast together with shapes (278,1) (280,2)" for FirstDifferenceOLS. The weights have to be finite as well as strictly positive, and the models raise "weights must be strictly positive and finite." for zero, negative and infinite weights. Missing weights are dropped with their observations, as before. The tests check the error of every model for -inf, inf, 0 and -1, and that missing weights give the model of the data without those observations. --- linearmodels/panel/model.py | 5 +- .../tests/panel/test_weights_validation.py | 80 +++++++++++++++++++ 2 files changed, 83 insertions(+), 2 deletions(-) create mode 100644 linearmodels/tests/panel/test_weights_validation.py diff --git a/linearmodels/panel/model.py b/linearmodels/panel/model.py index 9b635f8a7f..102b16cc10 100644 --- a/linearmodels/panel/model.py +++ b/linearmodels/panel/model.py @@ -575,8 +575,9 @@ def _validate_data(self) -> None: self._not_null = cast("BoolArray", np.asarray(~missing)) w_df = self.weights.dataframe - if np.any(np.asarray(w_df) <= 0): - raise ValueError("weights must be strictly positive.") + w_values = np.asarray(w_df) + if np.any(w_values <= 0) or np.any(np.isinf(w_values)): + raise ValueError("weights must be strictly positive and finite.") w_df = w_df / w_df.mean() self.weights = PanelData(w_df) rank_of_x = self._check_exog_rank() diff --git a/linearmodels/tests/panel/test_weights_validation.py b/linearmodels/tests/panel/test_weights_validation.py new file mode 100644 index 0000000000..e3e4c74eb4 --- /dev/null +++ b/linearmodels/tests/panel/test_weights_validation.py @@ -0,0 +1,80 @@ +""" +The weights of the panel models must be finite and strictly positive + +The weights are scaled to have a mean of one. Zero and negative weights were checked, but +an infinite weight was scaled to NaN and the other weights to zero, and the models failed +with errors from NumPy, SciPy and pandas that do not say what is wrong ("array must not +contain infs or NaNs", "SVD did not converge", an AbsorbingEffectError, and "operands could not +be broadcast together" for FirstDifferenceOLS). Missing weights are dropped with their +observations, as the other missing values are, and have to give the model of the data +without them. +""" + +import numpy as np +from numpy.testing import assert_allclose +import pandas as pd +import pytest + +from linearmodels.panel import ( + BetweenOLS, + FamaMacBeth, + FirstDifferenceOLS, + PanelOLS, + PooledOLS, + RandomEffects, +) +from linearmodels.shared.exceptions import MissingValueWarning + +MODELS = { + "pooled": lambda y, x, w: PooledOLS(y, x, weights=w), + "panel-entity-effects": lambda y, x, w: PanelOLS( + y, x, weights=w, entity_effects=True + ), + "panel-time-effects": lambda y, x, w: PanelOLS(y, x, weights=w, time_effects=True), + "between": lambda y, x, w: BetweenOLS(y, x, weights=w), + "first-difference": lambda y, x, w: FirstDifferenceOLS(y, x, weights=w), + "random-effects": lambda y, x, w: RandomEffects(y, x, weights=w), + "fama-macbeth": lambda y, x, w: FamaMacBeth(y, x, weights=w), +} + + +def make_panel(seed=0): + rg = np.random.default_rng(seed) + index = pd.MultiIndex.from_product( + [np.arange(30), np.arange(6)], names=["entity", "time"] + ) + x = pd.DataFrame( + rg.standard_normal((index.shape[0], 2)), index=index, columns=["a", "b"] + ) + y = pd.Series( + x @ [0.5, -0.3] + rg.standard_normal(index.shape[0]), index=index, name="y" + ) + w = pd.Series(rg.chisquare(4, index.shape[0]) / 4 + 0.2, index=index, name="w") + return y, x, w + + +@pytest.mark.parametrize("bad", [np.inf, -np.inf, 0.0, -1.0]) +@pytest.mark.parametrize("model", list(MODELS)) +def test_weights_must_be_finite_and_strictly_positive(model, bad): + y, x, w = make_panel() + w.iloc[7] = bad + with pytest.raises( + ValueError, match="weights must be strictly positive and finite" + ): + MODELS[model](y, x, w) + + +@pytest.mark.parametrize("model", list(MODELS)) +def test_missing_weights_are_dropped_with_their_observations(model): + y, x, w = make_panel() + missing = [3, 17, 18, 100] + w_missing = w.copy() + w_missing.iloc[missing] = np.nan + with pytest.warns(MissingValueWarning): + res = MODELS[model](y, x, w_missing).fit() + keep = np.ones(y.shape[0], dtype=bool) + keep[missing] = False + expected = MODELS[model](y[keep], x[keep], w[keep]).fit() + assert res.nobs == expected.nobs + assert_allclose(res.params, expected.params) + assert_allclose(res.std_errors, expected.std_errors)