flexcv.model_postprocessing

This module contains functions that are specifically written to fulfill the postprocessing requirements of the models. They are customized depending on what the model returns and what the user wants to log to Neptune. The functions are called by the cross_validate function in the flexcv.core module after the model has been fitted.

flexcv.model_postprocessing.LMERModelPostProcessor

Bases: ModelPostProcessor

Source code in flexcv/model_postprocessing.py

class LMERModelPostProcessor(ModelPostProcessor):
    def __init__(self):
        super().__init__()

    def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
        """Postprocessing function for the linear mixed effects model.
        Logs the summary of the model to neptune.

        Args:
            results_all_folds: A dict of results for all folds
            fold_result: A dataclass containing the results for the current fold
            run: neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            (dict): updated results dictionary
        """
        # LM is the only one where the regular logging wont work with get_params()
        # therefore, we need to get the summary via the model object and simply overwrite the empty logs
        params = fold_result.fit_result.get_summary()

        run[f"{fold_result.model_name}/Summary"].append(
            File.from_content(params, extension="html")
        )
        results_all_folds[fold_result.model_name]["parameters"][fold_result.k] = params
        return results_all_folds

flexcv.model_postprocessing.LMERModelPostProcessor.__call__(results_all_folds, fold_result, features, run, *args, **kwargs)

Postprocessing function for the linear mixed effects model. Logs the summary of the model to neptune.

Parameters:

Name	Type	Description	Default
results_all_folds		A dict of results for all folds	required
fold_result		A dataclass containing the results for the current fold	required
run		neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Type	Description
dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
    """Postprocessing function for the linear mixed effects model.
    Logs the summary of the model to neptune.

    Args:
        results_all_folds: A dict of results for all folds
        fold_result: A dataclass containing the results for the current fold
        run: neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        (dict): updated results dictionary
    """
    # LM is the only one where the regular logging wont work with get_params()
    # therefore, we need to get the summary via the model object and simply overwrite the empty logs
    params = fold_result.fit_result.get_summary()

    run[f"{fold_result.model_name}/Summary"].append(
        File.from_content(params, extension="html")
    )
    results_all_folds[fold_result.model_name]["parameters"][fold_result.k] = params
    return results_all_folds

flexcv.model_postprocessing.LinearModelPostProcessor

Bases: ModelPostProcessor

Source code in flexcv/model_postprocessing.py

class LinearModelPostProcessor(ModelPostProcessor):
    def __init__(self):
        super().__init__()

    def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
        """Postprocessing function for the linear regression model.
        Logs the summary of the model, the VIF and the plots to neptune.

        Args:
            results_all_folds: A dict of results for all folds
            fold_result: A dataclass containing the results for the current fold
            run: neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            (dict): updated results dictionary
        """
        # LM is the only one where the regular logging wont work with get_params()
        # therefore, we need to get the summary via the model object and simply overwrite the empty logs
        params = fold_result.fit_result.get_summary()

        run[f"{fold_result.model_name}/Summary"].append(
            File.from_content(params, extension="html")
        )
        results_all_folds[fold_result.model_name]["parameters"][fold_result.k] = params

        vif, fig, ax = LinearRegDiagnostic(fold_result.fit_result.md_)()  # type: ignore   instance has to be called after __init__
        run[f"{fold_result.model_name}/VIF"].append(File.as_html(vif))
        run[f"{fold_result.model_name}/Plots"].append(fig)
        del fig
        del vif
        del ax
        return results_all_folds

flexcv.model_postprocessing.LinearModelPostProcessor.__call__(results_all_folds, fold_result, features, run, *args, **kwargs)

Postprocessing function for the linear regression model. Logs the summary of the model, the VIF and the plots to neptune.

Parameters:

Name	Type	Description	Default
results_all_folds		A dict of results for all folds	required
fold_result		A dataclass containing the results for the current fold	required
run		neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Type	Description
dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
    """Postprocessing function for the linear regression model.
    Logs the summary of the model, the VIF and the plots to neptune.

    Args:
        results_all_folds: A dict of results for all folds
        fold_result: A dataclass containing the results for the current fold
        run: neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        (dict): updated results dictionary
    """
    # LM is the only one where the regular logging wont work with get_params()
    # therefore, we need to get the summary via the model object and simply overwrite the empty logs
    params = fold_result.fit_result.get_summary()

    run[f"{fold_result.model_name}/Summary"].append(
        File.from_content(params, extension="html")
    )
    results_all_folds[fold_result.model_name]["parameters"][fold_result.k] = params

    vif, fig, ax = LinearRegDiagnostic(fold_result.fit_result.md_)()  # type: ignore   instance has to be called after __init__
    run[f"{fold_result.model_name}/VIF"].append(File.as_html(vif))
    run[f"{fold_result.model_name}/Plots"].append(fig)
    del fig
    del vif
    del ax
    return results_all_folds

flexcv.model_postprocessing.LinearRegDiagnostic

Implementation by https://www.statsmodels.org/stable/examples/notebooks/generated/linear_regression_diagnostics_plots.html#Simple-multiple-linear-regression

Diagnostic plots to identify potential problems in a linear regression fit.

Mainly,

- non-linearity of data
- Correlation of error terms
- non-constant variance
- outliers
- high-leverage points
- collinearity

Authors

Prajwal Kafle (p33ajkafle@gmail.com, where 3 = r) Does not come with any sort of warranty. Please test the code one your end before using.

Matt Spinelli (m3spinelli@gmail.com, where 3 = r) (1) Fixed incorrect annotation of the top most extreme residuals in the Residuals vs Fitted and, especially, the Normal Q-Q plots. (2) Changed Residuals vs Leverage plot to match closer the y-axis range shown in the equivalent plot in the R package ggfortify. (3) Added horizontal line at y=0 in Residuals vs Leverage plot to match the plots in R package ggfortify and base R. (4) Added option for placing a vertical guideline on the Residuals vs Leverage plot using the rule of thumb of h = 2p/n to denote high leverage (high_leverage_threshold=True). (5) Added two more ways to compute the Cook's Distance (D) threshold: * 'baseR': D > 1 and D > 0.5 (default) * 'convention': D > 4/n * 'dof': D > 4 / (n - k - 1) (6) Fixed class name to conform to Pascal casing convention (7) Fixed Residuals vs Leverage legend to work with loc='best'

Source code in flexcv/model_postprocessing.py

class LinearRegDiagnostic:
    """
    Implementation by https://www.statsmodels.org/stable/examples/notebooks/generated/linear_regression_diagnostics_plots.html#Simple-multiple-linear-regression

    Diagnostic plots to identify potential problems in a linear regression fit.

    Mainly,

        - non-linearity of data
        - Correlation of error terms
        - non-constant variance
        - outliers
        - high-leverage points
        - collinearity

    Authors:
        Prajwal Kafle (p33ajkafle@gmail.com, where 3 = r)
        Does not come with any sort of warranty.
        Please test the code one your end before using.

        Matt Spinelli (m3spinelli@gmail.com, where 3 = r)
        (1) Fixed incorrect annotation of the top most extreme residuals in
            the Residuals vs Fitted and, especially, the Normal Q-Q plots.
        (2) Changed Residuals vs Leverage plot to match closer the y-axis
            range shown in the equivalent plot in the R package ggfortify.
        (3) Added horizontal line at y=0 in Residuals vs Leverage plot to
            match the plots in R package ggfortify and base R.
        (4) Added option for placing a vertical guideline on the Residuals
            vs Leverage plot using the rule of thumb of h = 2p/n to denote
            high leverage (high_leverage_threshold=True).
        (5) Added two more ways to compute the Cook's Distance (D) threshold:
            * 'baseR': D > 1 and D > 0.5 (default)
            * 'convention': D > 4/n
            * 'dof': D > 4 / (n - k - 1)
        (6) Fixed class name to conform to Pascal casing convention
        (7) Fixed Residuals vs Leverage legend to work with loc='best'
    """

    def __init__(
        self,
        results: Type[statsmodels.regression.linear_model.RegressionResultsWrapper],
    ) -> None:
        """
        For a linear regression model, generates following diagnostic plots:

        a. residual
        b. qq
        c. scale location and
        d. leverage

        and a table

        e. vif

        Args:
            results (Type[statsmodels.regression.linear_model.RegressionResultsWrapper]):
                must be instance of statsmodels.regression.linear_model object

        Raises:
            TypeError: if instance does not belong to above object

        Example:
        >>> import numpy as np
        >>> import pandas as pd
        >>> import statsmodels.formula.api as smf
        >>> x = np.linspace(-np.pi, np.pi, 100)
        >>> y = 3*x + 8 + np.random.normal(0,1, 100)
        >>> df = pd.DataFrame({'x':x, 'y':y})
        >>> res = smf.ols(formula= "y ~ x", data=df).fit()
        >>> cls = Linear_Reg_Diagnostic(res)
        >>> cls(plot_context="seaborn-paper")

        In case you do not need all plots you can also independently make an individual plot/table
        in following ways

        >>> cls = Linear_Reg_Diagnostic(res)
        >>> cls.residual_plot()
        >>> cls.qq_plot()
        >>> cls.scale_location_plot()
        >>> cls.leverage_plot()
        >>> cls.vif_table()
        """

        if (
            isinstance(
                results, statsmodels.regression.linear_model.RegressionResultsWrapper
            )
            is False
        ):
            raise TypeError(
                "result must be instance of statsmodels.regression.linear_model.RegressionResultsWrapper object"
            )

        self.results = maybe_unwrap_results(results)

        self.y_true = self.results.model.endog
        self.y_predict = self.results.fittedvalues
        self.xvar = self.results.model.exog
        self.xvar_names = self.results.model.exog_names

        self.residual = np.array(self.results.resid)
        influence = self.results.get_influence()
        self.residual_norm = influence.resid_studentized_internal
        self.leverage = influence.hat_matrix_diag
        self.cooks_distance = influence.cooks_distance[0]
        self.nparams = len(self.results.params)
        self.nresids = len(self.residual_norm)

    def __call__(self, plot_context="seaborn-v0_8-paper", **kwargs):
        # print(plt.style.available)
        with plt.style.context(plot_context):
            fig, ax = plt.subplots(nrows=2, ncols=2, figsize=(10, 10))
            self.residual_plot(ax=ax[0, 0])
            self.qq_plot(ax=ax[0, 1])
            self.scale_location_plot(ax=ax[1, 0])
            self.leverage_plot(
                ax=ax[1, 1],
                high_leverage_threshold=kwargs.get("high_leverage_threshold"),
                cooks_threshold=kwargs.get("cooks_threshold"),
            )
            # plt.show()

        return (
            self.vif_table(),
            fig,
            ax,
        )

    def residual_plot(self, ax=None):
        """
        Residual vs Fitted Plot

        Graphical tool to identify non-linearity.
        (Roughly) Horizontal red line is an indicator that the residual has a linear pattern
        """
        if ax is None:
            fig, ax = plt.subplots()

        sns.residplot(
            x=self.y_predict,
            y=self.residual,
            lowess=True,
            scatter_kws={"alpha": 0.5},
            line_kws={"color": "red", "lw": 1, "alpha": 0.8},
            ax=ax,
        )

        # annotations
        residual_abs = np.abs(self.residual)
        abs_resid = np.flip(np.argsort(residual_abs), 0)
        abs_resid_top_3 = abs_resid[:3]
        for i in abs_resid_top_3:
            ax.annotate(i, xy=(self.y_predict[i], self.residual[i]), color="C3")

        ax.set_title("Residuals vs Fitted", fontweight="bold")
        ax.set_xlabel("Fitted values")
        ax.set_ylabel("Residuals")
        return ax

    def qq_plot(self, ax=None):
        """
        Standarized Residual vs Theoretical Quantile plot

        Used to visually check if residuals are normally distributed.
        Points spread along the diagonal line will suggest so.
        """
        if ax is None:
            fig, ax = plt.subplots()

        QQ = ProbPlot(self.residual_norm)
        fig = QQ.qqplot(line="45", alpha=0.5, lw=1, ax=ax)

        # annotations
        abs_norm_resid = np.flip(np.argsort(np.abs(self.residual_norm)), 0)
        abs_norm_resid_top_3 = abs_norm_resid[:3]
        for i, x, y in self.__qq_top_resid(
            QQ.theoretical_quantiles, abs_norm_resid_top_3
        ):
            ax.annotate(i, xy=(x, y), ha="right", color="C3")

        ax.set_title("Normal Q-Q", fontweight="bold")
        ax.set_xlabel("Theoretical Quantiles")
        ax.set_ylabel("Standardized Residuals")
        return ax

    def scale_location_plot(self, ax=None):
        """
        Sqrt(Standarized Residual) vs Fitted values plot

        Used to check homoscedasticity of the residuals.
        Horizontal line will suggest so.
        """
        if ax is None:
            fig, ax = plt.subplots()

        residual_norm_abs_sqrt = np.sqrt(np.abs(self.residual_norm))

        ax.scatter(self.y_predict, residual_norm_abs_sqrt, alpha=0.5)
        sns.regplot(
            x=self.y_predict,
            y=residual_norm_abs_sqrt,
            scatter=False,
            ci=False,
            lowess=True,
            line_kws={"color": "red", "lw": 1, "alpha": 0.8},
            ax=ax,
        )

        # annotations
        abs_sq_norm_resid = np.flip(np.argsort(residual_norm_abs_sqrt), 0)
        abs_sq_norm_resid_top_3 = abs_sq_norm_resid[:3]
        for i in abs_sq_norm_resid_top_3:
            ax.annotate(
                i, xy=(self.y_predict[i], residual_norm_abs_sqrt[i]), color="C3"
            )

        ax.set_title("Scale-Location", fontweight="bold")
        ax.set_xlabel("Fitted values")
        ax.set_ylabel(r"$\sqrt{|\mathrm{Standardized\ Residuals}|}$")
        return ax

    def leverage_plot(
        self, ax=None, high_leverage_threshold=False, cooks_threshold="baseR"
    ):
        """
        Residual vs Leverage plot

        Points falling outside Cook's distance curves are considered observation that can sway the fit
        aka are influential.
        Good to have none outside the curves.
        """
        if ax is None:
            fig, ax = plt.subplots()

        ax.scatter(self.leverage, self.residual_norm, alpha=0.5)

        sns.regplot(
            x=self.leverage,
            y=self.residual_norm,
            scatter=False,
            ci=False,
            lowess=True,
            line_kws={"color": "red", "lw": 1, "alpha": 0.8},
            ax=ax,
        )

        # annotations
        leverage_top_3 = np.flip(np.argsort(self.cooks_distance), 0)[:3]
        for i in leverage_top_3:
            ax.annotate(i, xy=(self.leverage[i], self.residual_norm[i]), color="C3")

        factors = []
        if cooks_threshold == "baseR" or cooks_threshold is None:
            factors = [1, 0.5]
        elif cooks_threshold == "convention":
            factors = [4 / self.nresids]
        elif cooks_threshold == "dof":
            factors = [4 / (self.nresids - self.nparams)]
        else:
            raise ValueError(
                "threshold_method must be one of the following: 'convention', 'dof', or 'baseR' (default)"
            )
        for i, factor in enumerate(factors):
            label = "Cook's distance" if i == 0 else None
            xtemp, ytemp = self.__cooks_dist_line(factor)
            ax.plot(xtemp, ytemp, label=label, lw=1.25, ls="--", color="red")
            ax.plot(xtemp, np.negative(ytemp), lw=1.25, ls="--", color="red")

        if high_leverage_threshold:
            high_leverage = 2 * self.nparams / self.nresids
            if max(self.leverage) > high_leverage:
                ax.axvline(
                    high_leverage, label="High leverage", ls="-.", color="purple", lw=1
                )

        ax.axhline(0, ls="dotted", color="black", lw=1.25)
        ax.set_xlim(0, max(self.leverage) + 0.01)
        ax.set_ylim(min(self.residual_norm) - 0.1, max(self.residual_norm) + 0.1)
        ax.set_title("Residuals vs Leverage", fontweight="bold")
        ax.set_xlabel("Leverage")
        ax.set_ylabel("Standardized Residuals")
        plt.legend(loc="best")
        return ax

    def vif_table(self):
        """
        VIF table

        VIF, the variance inflation factor, is a measure of multicollinearity.
        VIF > 5 for a variable indicates that it is highly collinear with the
        other input variables.
        """
        vif_df = pd.DataFrame()
        vif_df["Features"] = self.xvar_names
        vif_df["VIF Factor"] = [
            variance_inflation_factor(self.xvar, i) for i in range(self.xvar.shape[1])
        ]

        return vif_df.sort_values("VIF Factor").round(2)

    def __cooks_dist_line(self, factor):
        """
        Helper function for plotting Cook's distance curves
        """
        p = self.nparams
        formula = lambda x: np.sqrt((factor * p * (1 - x)) / x)
        x = np.linspace(0.001, max(self.leverage), 50)
        y = formula(x)
        return x, y

    def __qq_top_resid(self, quantiles, top_residual_indices):
        """
        Helper generator function yielding the index and coordinates
        """
        offset = 0
        quant_index = 0
        previous_is_negative = None
        for resid_index in top_residual_indices:
            y = self.residual_norm[resid_index]
            is_negative = y < 0
            if previous_is_negative == None or previous_is_negative == is_negative:
                offset += 1
            else:
                quant_index -= offset
            x = (
                quantiles[quant_index]
                if is_negative
                else np.flip(quantiles, 0)[quant_index]
            )
            quant_index += 1
            previous_is_negative = is_negative
            yield resid_index, x, y

flexcv.model_postprocessing.LinearRegDiagnostic.__cooks_dist_line(factor)

Helper function for plotting Cook's distance curves

Source code in flexcv/model_postprocessing.py

def __cooks_dist_line(self, factor):
    """
    Helper function for plotting Cook's distance curves
    """
    p = self.nparams
    formula = lambda x: np.sqrt((factor * p * (1 - x)) / x)
    x = np.linspace(0.001, max(self.leverage), 50)
    y = formula(x)
    return x, y

flexcv.model_postprocessing.LinearRegDiagnostic.__init__(results)

For a linear regression model, generates following diagnostic plots:

a. residual b. qq c. scale location and d. leverage

and a table

e. vif

Parameters:

Name	Type	Description	Default
results	Type[RegressionResultsWrapper]	must be instance of statsmodels.regression.linear_model object	required

Raises:

Type	Description
TypeError	if instance does not belong to above object

Example:

import numpy as np import pandas as pd import statsmodels.formula.api as smf x = np.linspace(-np.pi, np.pi, 100) y = 3*x + 8 + np.random.normal(0,1, 100) df = pd.DataFrame({'x':x, 'y':y}) res = smf.ols(formula= "y ~ x", data=df).fit() cls = Linear_Reg_Diagnostic(res) cls(plot_context="seaborn-paper")

In case you do not need all plots you can also independently make an individual plot/table in following ways

cls = Linear_Reg_Diagnostic(res) cls.residual_plot() cls.qq_plot() cls.scale_location_plot() cls.leverage_plot() cls.vif_table()

Source code in flexcv/model_postprocessing.py

def __init__(
    self,
    results: Type[statsmodels.regression.linear_model.RegressionResultsWrapper],
) -> None:
    """
    For a linear regression model, generates following diagnostic plots:

    a. residual
    b. qq
    c. scale location and
    d. leverage

    and a table

    e. vif

    Args:
        results (Type[statsmodels.regression.linear_model.RegressionResultsWrapper]):
            must be instance of statsmodels.regression.linear_model object

    Raises:
        TypeError: if instance does not belong to above object

    Example:
    >>> import numpy as np
    >>> import pandas as pd
    >>> import statsmodels.formula.api as smf
    >>> x = np.linspace(-np.pi, np.pi, 100)
    >>> y = 3*x + 8 + np.random.normal(0,1, 100)
    >>> df = pd.DataFrame({'x':x, 'y':y})
    >>> res = smf.ols(formula= "y ~ x", data=df).fit()
    >>> cls = Linear_Reg_Diagnostic(res)
    >>> cls(plot_context="seaborn-paper")

    In case you do not need all plots you can also independently make an individual plot/table
    in following ways

    >>> cls = Linear_Reg_Diagnostic(res)
    >>> cls.residual_plot()
    >>> cls.qq_plot()
    >>> cls.scale_location_plot()
    >>> cls.leverage_plot()
    >>> cls.vif_table()
    """

    if (
        isinstance(
            results, statsmodels.regression.linear_model.RegressionResultsWrapper
        )
        is False
    ):
        raise TypeError(
            "result must be instance of statsmodels.regression.linear_model.RegressionResultsWrapper object"
        )

    self.results = maybe_unwrap_results(results)

    self.y_true = self.results.model.endog
    self.y_predict = self.results.fittedvalues
    self.xvar = self.results.model.exog
    self.xvar_names = self.results.model.exog_names

    self.residual = np.array(self.results.resid)
    influence = self.results.get_influence()
    self.residual_norm = influence.resid_studentized_internal
    self.leverage = influence.hat_matrix_diag
    self.cooks_distance = influence.cooks_distance[0]
    self.nparams = len(self.results.params)
    self.nresids = len(self.residual_norm)

flexcv.model_postprocessing.LinearRegDiagnostic.__qq_top_resid(quantiles, top_residual_indices)

Helper generator function yielding the index and coordinates

Source code in flexcv/model_postprocessing.py

def __qq_top_resid(self, quantiles, top_residual_indices):
    """
    Helper generator function yielding the index and coordinates
    """
    offset = 0
    quant_index = 0
    previous_is_negative = None
    for resid_index in top_residual_indices:
        y = self.residual_norm[resid_index]
        is_negative = y < 0
        if previous_is_negative == None or previous_is_negative == is_negative:
            offset += 1
        else:
            quant_index -= offset
        x = (
            quantiles[quant_index]
            if is_negative
            else np.flip(quantiles, 0)[quant_index]
        )
        quant_index += 1
        previous_is_negative = is_negative
        yield resid_index, x, y

flexcv.model_postprocessing.LinearRegDiagnostic.leverage_plot(ax=None, high_leverage_threshold=False, cooks_threshold='baseR')

Residual vs Leverage plot

Points falling outside Cook's distance curves are considered observation that can sway the fit aka are influential. Good to have none outside the curves.

Source code in flexcv/model_postprocessing.py

def leverage_plot(
    self, ax=None, high_leverage_threshold=False, cooks_threshold="baseR"
):
    """
    Residual vs Leverage plot

    Points falling outside Cook's distance curves are considered observation that can sway the fit
    aka are influential.
    Good to have none outside the curves.
    """
    if ax is None:
        fig, ax = plt.subplots()

    ax.scatter(self.leverage, self.residual_norm, alpha=0.5)

    sns.regplot(
        x=self.leverage,
        y=self.residual_norm,
        scatter=False,
        ci=False,
        lowess=True,
        line_kws={"color": "red", "lw": 1, "alpha": 0.8},
        ax=ax,
    )

    # annotations
    leverage_top_3 = np.flip(np.argsort(self.cooks_distance), 0)[:3]
    for i in leverage_top_3:
        ax.annotate(i, xy=(self.leverage[i], self.residual_norm[i]), color="C3")

    factors = []
    if cooks_threshold == "baseR" or cooks_threshold is None:
        factors = [1, 0.5]
    elif cooks_threshold == "convention":
        factors = [4 / self.nresids]
    elif cooks_threshold == "dof":
        factors = [4 / (self.nresids - self.nparams)]
    else:
        raise ValueError(
            "threshold_method must be one of the following: 'convention', 'dof', or 'baseR' (default)"
        )
    for i, factor in enumerate(factors):
        label = "Cook's distance" if i == 0 else None
        xtemp, ytemp = self.__cooks_dist_line(factor)
        ax.plot(xtemp, ytemp, label=label, lw=1.25, ls="--", color="red")
        ax.plot(xtemp, np.negative(ytemp), lw=1.25, ls="--", color="red")

    if high_leverage_threshold:
        high_leverage = 2 * self.nparams / self.nresids
        if max(self.leverage) > high_leverage:
            ax.axvline(
                high_leverage, label="High leverage", ls="-.", color="purple", lw=1
            )

    ax.axhline(0, ls="dotted", color="black", lw=1.25)
    ax.set_xlim(0, max(self.leverage) + 0.01)
    ax.set_ylim(min(self.residual_norm) - 0.1, max(self.residual_norm) + 0.1)
    ax.set_title("Residuals vs Leverage", fontweight="bold")
    ax.set_xlabel("Leverage")
    ax.set_ylabel("Standardized Residuals")
    plt.legend(loc="best")
    return ax

flexcv.model_postprocessing.LinearRegDiagnostic.qq_plot(ax=None)

Standarized Residual vs Theoretical Quantile plot

Used to visually check if residuals are normally distributed. Points spread along the diagonal line will suggest so.

Source code in flexcv/model_postprocessing.py

def qq_plot(self, ax=None):
    """
    Standarized Residual vs Theoretical Quantile plot

    Used to visually check if residuals are normally distributed.
    Points spread along the diagonal line will suggest so.
    """
    if ax is None:
        fig, ax = plt.subplots()

    QQ = ProbPlot(self.residual_norm)
    fig = QQ.qqplot(line="45", alpha=0.5, lw=1, ax=ax)

    # annotations
    abs_norm_resid = np.flip(np.argsort(np.abs(self.residual_norm)), 0)
    abs_norm_resid_top_3 = abs_norm_resid[:3]
    for i, x, y in self.__qq_top_resid(
        QQ.theoretical_quantiles, abs_norm_resid_top_3
    ):
        ax.annotate(i, xy=(x, y), ha="right", color="C3")

    ax.set_title("Normal Q-Q", fontweight="bold")
    ax.set_xlabel("Theoretical Quantiles")
    ax.set_ylabel("Standardized Residuals")
    return ax

flexcv.model_postprocessing.LinearRegDiagnostic.residual_plot(ax=None)

Residual vs Fitted Plot

Graphical tool to identify non-linearity. (Roughly) Horizontal red line is an indicator that the residual has a linear pattern

Source code in flexcv/model_postprocessing.py

def residual_plot(self, ax=None):
    """
    Residual vs Fitted Plot

    Graphical tool to identify non-linearity.
    (Roughly) Horizontal red line is an indicator that the residual has a linear pattern
    """
    if ax is None:
        fig, ax = plt.subplots()

    sns.residplot(
        x=self.y_predict,
        y=self.residual,
        lowess=True,
        scatter_kws={"alpha": 0.5},
        line_kws={"color": "red", "lw": 1, "alpha": 0.8},
        ax=ax,
    )

    # annotations
    residual_abs = np.abs(self.residual)
    abs_resid = np.flip(np.argsort(residual_abs), 0)
    abs_resid_top_3 = abs_resid[:3]
    for i in abs_resid_top_3:
        ax.annotate(i, xy=(self.y_predict[i], self.residual[i]), color="C3")

    ax.set_title("Residuals vs Fitted", fontweight="bold")
    ax.set_xlabel("Fitted values")
    ax.set_ylabel("Residuals")
    return ax

flexcv.model_postprocessing.LinearRegDiagnostic.scale_location_plot(ax=None)

Sqrt(Standarized Residual) vs Fitted values plot

Used to check homoscedasticity of the residuals. Horizontal line will suggest so.

Source code in flexcv/model_postprocessing.py

def scale_location_plot(self, ax=None):
    """
    Sqrt(Standarized Residual) vs Fitted values plot

    Used to check homoscedasticity of the residuals.
    Horizontal line will suggest so.
    """
    if ax is None:
        fig, ax = plt.subplots()

    residual_norm_abs_sqrt = np.sqrt(np.abs(self.residual_norm))

    ax.scatter(self.y_predict, residual_norm_abs_sqrt, alpha=0.5)
    sns.regplot(
        x=self.y_predict,
        y=residual_norm_abs_sqrt,
        scatter=False,
        ci=False,
        lowess=True,
        line_kws={"color": "red", "lw": 1, "alpha": 0.8},
        ax=ax,
    )

    # annotations
    abs_sq_norm_resid = np.flip(np.argsort(residual_norm_abs_sqrt), 0)
    abs_sq_norm_resid_top_3 = abs_sq_norm_resid[:3]
    for i in abs_sq_norm_resid_top_3:
        ax.annotate(
            i, xy=(self.y_predict[i], residual_norm_abs_sqrt[i]), color="C3"
        )

    ax.set_title("Scale-Location", fontweight="bold")
    ax.set_xlabel("Fitted values")
    ax.set_ylabel(r"$\sqrt{|\mathrm{Standardized\ Residuals}|}$")
    return ax

flexcv.model_postprocessing.LinearRegDiagnostic.vif_table()

VIF table

VIF, the variance inflation factor, is a measure of multicollinearity. VIF > 5 for a variable indicates that it is highly collinear with the other input variables.

Source code in flexcv/model_postprocessing.py

def vif_table(self):
    """
    VIF table

    VIF, the variance inflation factor, is a measure of multicollinearity.
    VIF > 5 for a variable indicates that it is highly collinear with the
    other input variables.
    """
    vif_df = pd.DataFrame()
    vif_df["Features"] = self.xvar_names
    vif_df["VIF Factor"] = [
        variance_inflation_factor(self.xvar, i) for i in range(self.xvar.shape[1])
    ]

    return vif_df.sort_values("VIF Factor").round(2)

flexcv.model_postprocessing.MERFModelPostProcessor

Bases: ModelPostProcessor

Source code in flexcv/model_postprocessing.py

class MERFModelPostProcessor(ModelPostProcessor):
    def __init__(self):
        super().__init__()

    def __call__(
        self, results_all_folds, fold_result, y_pred_base, run, *args, **kwargs
    ):
        """Postprocessing function for the expectation maximization model (MERF).
        Logs training and test plots to Neptune.

        Args:
            results_all_folds: A dict of results for all folds
            fold_result: A dataclass containing the results for the current fold
            run: neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            (dict): updated results dictionary
        """
        run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
            {
                "max_iterations": fold_result.fit_result.max_iterations,
            }
        )

        plot.plot_merf_results(
            y=fold_result.y_test,
            yhat_base=y_pred_base,
            yhat_me_model=fold_result.y_pred,
            model=fold_result.fit_result,
            model_name=kwargs["mixed_name"],
            max_iterations=fold_result.fit_result.max_iterations,
            run=run,
        )
        plt.close()
        return results_all_folds

flexcv.model_postprocessing.MERFModelPostProcessor.__call__(results_all_folds, fold_result, y_pred_base, run, *args, **kwargs)

Postprocessing function for the expectation maximization model (MERF). Logs training and test plots to Neptune.

Parameters:

Name	Type	Description	Default
results_all_folds		A dict of results for all folds	required
fold_result		A dataclass containing the results for the current fold	required
run		neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Type	Description
dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(
    self, results_all_folds, fold_result, y_pred_base, run, *args, **kwargs
):
    """Postprocessing function for the expectation maximization model (MERF).
    Logs training and test plots to Neptune.

    Args:
        results_all_folds: A dict of results for all folds
        fold_result: A dataclass containing the results for the current fold
        run: neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        (dict): updated results dictionary
    """
    run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
        {
            "max_iterations": fold_result.fit_result.max_iterations,
        }
    )

    plot.plot_merf_results(
        y=fold_result.y_test,
        yhat_base=y_pred_base,
        yhat_me_model=fold_result.y_pred,
        model=fold_result.fit_result,
        model_name=kwargs["mixed_name"],
        max_iterations=fold_result.fit_result.max_iterations,
        run=run,
    )
    plt.close()
    return results_all_folds

flexcv.model_postprocessing.ModelPostProcessor

Bases: ABC

Abstract base class for model postprocessing functions. All postprocessing functions must inherit from this class. Implement your own post processing routing by inheriting from this class and implementing the call method. The class instance is called in the cross validation loop.

Methods:

Name	Description
__call__	method to be implemented by the user. This method is called in the cross validation loop.

Source code in flexcv/model_postprocessing.py

class ModelPostProcessor(ABC):
    """Abstract base class for model postprocessing functions.
    All postprocessing functions must inherit from this class.
    Implement your own post processing routing by inheriting from this class and implementing the __call__ method.
    The class instance is called in the cross validation loop.

    Methods:
        __call__ (abstract): method to be implemented by the user. This method is called in the cross validation loop.
    """

    def __init__(self):
        pass

    def __call__(
        self,
        results_all_folds: dict,
        fold_result: SingleModelFoldResult,
        features: pd.Index | list[str] | np.ndarray[str],
        run: Run,
        *args,
        **kwargs,
    ) -> dict:
        """This method is called in the cross validation loop.
        Implement your own post processing routing by inheriting from this class and implementing the __call__ method.

        Args:
            results_all_folds (dict): A dict of results for all folds
            fold_result (SingleModelFoldResult): A dataclass containing the results for the current fold
            features (pd.Index | list[str] | np.ndarray[str]): The features used in the model
            run (Run): neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            results_all_folds (dict): updated results dictionary
        """

        pass

flexcv.model_postprocessing.ModelPostProcessor.__call__(results_all_folds, fold_result, features, run, *args, **kwargs)

This method is called in the cross validation loop. Implement your own post processing routing by inheriting from this class and implementing the call method.

Parameters:

Name	Type	Description	Default
results_all_folds	dict	A dict of results for all folds	required
fold_result	SingleModelFoldResult	A dataclass containing the results for the current fold	required
features	Index | list[str] | ndarray[str]	The features used in the model	required
run	Run	neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Name	Type	Description
results_all_folds	dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(
    self,
    results_all_folds: dict,
    fold_result: SingleModelFoldResult,
    features: pd.Index | list[str] | np.ndarray[str],
    run: Run,
    *args,
    **kwargs,
) -> dict:
    """This method is called in the cross validation loop.
    Implement your own post processing routing by inheriting from this class and implementing the __call__ method.

    Args:
        results_all_folds (dict): A dict of results for all folds
        fold_result (SingleModelFoldResult): A dataclass containing the results for the current fold
        features (pd.Index | list[str] | np.ndarray[str]): The features used in the model
        run (Run): neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        results_all_folds (dict): updated results dictionary
    """

    pass

flexcv.model_postprocessing.RandomForestModelPostProcessor

Bases: ModelPostProcessor

Source code in flexcv/model_postprocessing.py

class RandomForestModelPostProcessor(ModelPostProcessor):
    def __init__(self):
        super().__init__()

    def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
        """Postprocessing function for the random forest model.
        Logs the parameters to Neptune.
        Generates beeswarm plots from SHAP explainers for the training and test data and logs them to Neptune.

        Args:
            results_all_folds: A dict of results for all folds
            fold_result: A dataclass containing the results for the current fold
            run: neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            (dict): updated results dictionary
        """
        run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
            fold_result.best_params
        )
        explainer = shap.TreeExplainer(fold_result.best_model)
        shap_values = explainer.shap_values(fold_result.X_train)
        plot.plot_shap(
            shap_values=shap_values,
            X=fold_result.X_train,
            run=run,
            log_destination=f"{fold_result.model_name}/SHAP/Fold",
            dependency=False,
        )
        if "shap_values" not in results_all_folds[fold_result.model_name].keys():
            results_all_folds[fold_result.model_name]["shap_values"] = [shap_values]
        else:
            results_all_folds[fold_result.model_name]["shap_values"].append(shap_values)

        shap_values_test = explainer.shap_values(fold_result.X_test)
        plot.plot_shap(
            shap_values=shap_values_test,
            X=fold_result.X_test,
            run=run,
            log_destination=f"{fold_result.model_name}/SHAP/Test_Fold",
            dependency=False,
        )

        return results_all_folds

flexcv.model_postprocessing.RandomForestModelPostProcessor.__call__(results_all_folds, fold_result, features, run, *args, **kwargs)

Postprocessing function for the random forest model. Logs the parameters to Neptune. Generates beeswarm plots from SHAP explainers for the training and test data and logs them to Neptune.

Parameters:

Name	Type	Description	Default
results_all_folds		A dict of results for all folds	required
fold_result		A dataclass containing the results for the current fold	required
run		neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Type	Description
dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
    """Postprocessing function for the random forest model.
    Logs the parameters to Neptune.
    Generates beeswarm plots from SHAP explainers for the training and test data and logs them to Neptune.

    Args:
        results_all_folds: A dict of results for all folds
        fold_result: A dataclass containing the results for the current fold
        run: neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        (dict): updated results dictionary
    """
    run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
        fold_result.best_params
    )
    explainer = shap.TreeExplainer(fold_result.best_model)
    shap_values = explainer.shap_values(fold_result.X_train)
    plot.plot_shap(
        shap_values=shap_values,
        X=fold_result.X_train,
        run=run,
        log_destination=f"{fold_result.model_name}/SHAP/Fold",
        dependency=False,
    )
    if "shap_values" not in results_all_folds[fold_result.model_name].keys():
        results_all_folds[fold_result.model_name]["shap_values"] = [shap_values]
    else:
        results_all_folds[fold_result.model_name]["shap_values"].append(shap_values)

    shap_values_test = explainer.shap_values(fold_result.X_test)
    plot.plot_shap(
        shap_values=shap_values_test,
        X=fold_result.X_test,
        run=run,
        log_destination=f"{fold_result.model_name}/SHAP/Test_Fold",
        dependency=False,
    )

    return results_all_folds

flexcv.model_postprocessing.SVRModelPostProcessor

Bases: ModelPostProcessor

Source code in flexcv/model_postprocessing.py

class SVRModelPostProcessor(ModelPostProcessor):
    def __init__(self):
        super().__init__()

    def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
        """Postprocessing function for the SVR model.
        Logs the parameters to Neptune.
        Logs permutation importance plots to Neptune.

        Args:
            results_all_folds: A dict of results for all folds
            fold_result: A dataclass containing the results for the current fold
            run: neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            (dict): updated results dictionary
        """
        with plt.style.context("ggplot"):
            run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
                fold_result.best_params
            )

            fig, tmp_df = permutation_importance(
                fold_result.best_model,
                fold_result.model_name,
                fold_result.X_test,
                fold_result.y_test,
                fold_result.X_test.columns,
            )
            run[f"{fold_result.model_name}/PermFeatImportance/Figures"].append(fig)
            run[f"{fold_result.model_name}/PermFeatImportance/Table"].append(File.as_html(tmp_df))

        return results_all_folds

flexcv.model_postprocessing.SVRModelPostProcessor.__call__(results_all_folds, fold_result, features, run, *args, **kwargs)

Postprocessing function for the SVR model. Logs the parameters to Neptune. Logs permutation importance plots to Neptune.

Parameters:

Name	Type	Description	Default
results_all_folds		A dict of results for all folds	required
fold_result		A dataclass containing the results for the current fold	required
run		neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Type	Description
dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
    """Postprocessing function for the SVR model.
    Logs the parameters to Neptune.
    Logs permutation importance plots to Neptune.

    Args:
        results_all_folds: A dict of results for all folds
        fold_result: A dataclass containing the results for the current fold
        run: neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        (dict): updated results dictionary
    """
    with plt.style.context("ggplot"):
        run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
            fold_result.best_params
        )

        fig, tmp_df = permutation_importance(
            fold_result.best_model,
            fold_result.model_name,
            fold_result.X_test,
            fold_result.y_test,
            fold_result.X_test.columns,
        )
        run[f"{fold_result.model_name}/PermFeatImportance/Figures"].append(fig)
        run[f"{fold_result.model_name}/PermFeatImportance/Table"].append(File.as_html(tmp_df))

    return results_all_folds

flexcv.model_postprocessing.XGBoostModelPostProcessor

Bases: ModelPostProcessor

Source code in flexcv/model_postprocessing.py

class XGBoostModelPostProcessor(ModelPostProcessor):
    def __init__(self):
        super().__init__()

    def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
        """Postprocessing function for the xgboost model.
        Logs the parameters to Neptune.
        Generates beeswarm plots from SHAP explainers for the training and test data and logs them to Neptune.

        Args:
            results_all_folds: A dict of results for all folds
            fold_result: A dataclass containing the results for the current fold
            run: neptune run object
            *args: any additional arguments
            **kwargs: any additional keyword arguments

        Returns:
            (dict): updated results dictionary
        """
        explainer = shap.TreeExplainer(fold_result.best_model)
        shap_values = explainer.shap_values(fold_result.X_train)
        plot.plot_shap(
            shap_values=shap_values,
            X=fold_result.X_train,
            run=run,
            log_destination=f"{fold_result.model_name}/SHAP/Train_Fold",
            dependency=False,
        )
        if "shap_values" not in results_all_folds[fold_result.model_name].keys():
            results_all_folds[fold_result.model_name]["shap_values"] = [shap_values]
        else:
            results_all_folds[fold_result.model_name]["shap_values"].append(shap_values)

        run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
            fold_result.best_params
        )

        shap_values_test = explainer.shap_values(fold_result.X_test)
        plot.plot_shap(
            shap_values=shap_values_test,
            X=fold_result.X_test,
            run=run,
            log_destination=f"{fold_result.model_name}/SHAP/Test_Fold",
            dependency=False,
        )

        return results_all_folds

flexcv.model_postprocessing.XGBoostModelPostProcessor.__call__(results_all_folds, fold_result, features, run, *args, **kwargs)

Postprocessing function for the xgboost model. Logs the parameters to Neptune. Generates beeswarm plots from SHAP explainers for the training and test data and logs them to Neptune.

Parameters:

Name	Type	Description	Default
results_all_folds		A dict of results for all folds	required
fold_result		A dataclass containing the results for the current fold	required
run		neptune run object	required
*args		any additional arguments	()
**kwargs		any additional keyword arguments	{}

Returns:

Type	Description
dict	updated results dictionary

Source code in flexcv/model_postprocessing.py

def __call__(self, results_all_folds, fold_result, features, run, *args, **kwargs):
    """Postprocessing function for the xgboost model.
    Logs the parameters to Neptune.
    Generates beeswarm plots from SHAP explainers for the training and test data and logs them to Neptune.

    Args:
        results_all_folds: A dict of results for all folds
        fold_result: A dataclass containing the results for the current fold
        run: neptune run object
        *args: any additional arguments
        **kwargs: any additional keyword arguments

    Returns:
        (dict): updated results dictionary
    """
    explainer = shap.TreeExplainer(fold_result.best_model)
    shap_values = explainer.shap_values(fold_result.X_train)
    plot.plot_shap(
        shap_values=shap_values,
        X=fold_result.X_train,
        run=run,
        log_destination=f"{fold_result.model_name}/SHAP/Train_Fold",
        dependency=False,
    )
    if "shap_values" not in results_all_folds[fold_result.model_name].keys():
        results_all_folds[fold_result.model_name]["shap_values"] = [shap_values]
    else:
        results_all_folds[fold_result.model_name]["shap_values"].append(shap_values)

    run[f"{fold_result.model_name}/BestParams/{fold_result.k}"] = pformat(
        fold_result.best_params
    )

    shap_values_test = explainer.shap_values(fold_result.X_test)
    plot.plot_shap(
        shap_values=shap_values_test,
        X=fold_result.X_test,
        run=run,
        log_destination=f"{fold_result.model_name}/SHAP/Test_Fold",
        dependency=False,
    )

    return results_all_folds