Suggestion: New Constraint, Positive Only

[FloeFoxon]

I would like to share a suggestion for improved functionality, please. pyGAM seems to have almost all of the functionality I need for my particular use case, which is extraordinary, as comparator packages don't even come close. However, in my use case, I have contextual reasons why the spline term must always be positive (≥0).

As an example, please see the figure below, which shows the results of fitting a model with a single spline term and no intercept, and which was generated with the code

XX = gam.generate_X_grid(term=0)
ax.plot(XX[:, 0], gam.partial_dependence(term=0, X=XX))

similar to the code provided in the 'Tour of pyGAM' in the documentation.

The spline term is fitted with splines and is constrained to be monotonic and concave. The shape is correct. However, the splines values are positive for x-axis values greater than about 3, but are negative for x-axis values less than about 3. In my use case, it is necessary to further constrain the splines so that they cannot fall below zero, i.e. they must be positive only. Ideally, one could then fit a model with code for the spline term something like

s(0, n_splines = 10, constraints=['concave', 'monotonic_inc', 'positive_only'])

I have not yet delved into the source code to see how easy this would be to implement, but I intend to do so, unless the authors of pyGAM are able to implement this suggestion more quickly and painlessly than I could. Thanks so much for the package and for reading this suggestion.

[FloeFoxon]

Well, I took a look at the source code and it didn't take long to come up with the following. The required changes are,

(1) The addition of the following function to the penalties.py file:

def positive_only(n, coef):
    """
    Builds a penalty matrix for P-Splines with positive features.
    Penalizes violation of a positive feature function.

    Parameters
    ----------
    n : int
        number of splines
    coef : array-like, coefficients of the feature function

    Returns
    -------
    penalty matrix : sparse csc matrix of shape (n,n)
    """
    return sp.sparse.csc_matrix(np.diag((coef.ravel() < 0).astype(float)))

(2) Editing the CONSTRAINTS dictionary in the same penalties.py file to:

CONSTRAINTS = {
    'convex': convex,
    'concave': concave,
    'monotonic_inc': monotonic_inc,
    'monotonic_dec': monotonic_dec,
    'positive_only': positive_only,
    'none': none,
}

The code I developed above seems to work for my use case. This figure below shows the same model and data as in my previous comment using s(0, n_splines = 10, constraints=['concave', 'monotonic_inc', 'positive_only']). As you can see, it no longer falls below zero, as required: