# Linear Models

There's a few linear models out there that we felt were generally useful. This
document will highlight some of them. 

## LOWESS

Lowess stands for LOcally WEighted regreSSion and has historically 
been used for smoothing but you can also use it
for machine learning where you're interested in interpolating. Here's 
a demo; 


```python
import numpy as np
import matplotlib.pylab as plt
from sklego.linear_model import LowessRegression

n = 100
xs = np.linspace(0, np.pi, n)
ys = 1 + np.sin(xs) + np.cos(xs**2) + np.random.normal(0, 0.1, n)

mod = LowessRegression(sigma=0.1, span=0.5).fit(xs.reshape(-1, 1), ys)

xs_new = np.linspace(-1, np.pi + 1, n*2)
preds = mod.predict(xs_new.reshape(-1, 1))

plt.figure(figsize=(12, 4))
plt.scatter(xs, ys)
plt.plot(xs_new, preds, color='orange')
plt.title("Be careful with extrapolation here.");
```

The line does not look linear but that's because internally, during prediction,
many weighted linear regressions are happening. The gif below demonstrates how 
the data is being weighted when we would make a prediction. 



### Details on `sigma`.

We'll also show two different prediction outcomes depending on the hyperparameter
`sigma`. 




You may be tempted now to think that a lower sigma always has a better fit, but you
need to be careful here. The data might have gaps and larger sigma values will be
able to properly regularize.



Note that this regression also works in higher dimensions but the main downside of 
this approach is that it is really slow when making predictions. 

If you really want to get advanced there's also a hyperparameter `span` but you'll
really need to know what you're doing. It was added for completeness but the authors 
of this package have yet to find a proper usecase for it. 

### Details on `span`.

The `span` parameter can be used to force that you'll only include a certain
  percentage of your dataset. Technically without a `span` you are still using
  all the data in your dataset, albeit with small weights if they are far away.
  
The effect of the `span` parameter on the weights can be seen below.


  
This will also effect the predictions.</p>



You may need to squint your eyes a bit to see it, but lower spans cause more 
  jiggles and less smooth curves.

## ProbWeightRegression

> Note, this is a somewhat experimental feature. We found this feature to be plausibly useful but we've not seen it cause a big "win" yet. 

Let's say that you're interested in combining a few models for a regression task. 
You could put them together in an ensemble. Say we've got predictions $y_1, y_2, y_3$, each of which come from respectable models, then you may want to combine these together in another linear regression. This way the new, hopefully best, prediction $y_*$ is defined via;

$$y^* = w_1 y_1 + w_2 y_2 + w_3 y_3$$ 

This can be a valid way of reweighing. But there's one issue:
technically the weights $w_1, w_2, w_3$ can sum to a number that isn't one. Since that's 
numerically totally possible we need to be aware that we can end up in a strange situation. 

The `ProbWeightRegression` adresses this. It assumes that every input it receives is
the output of a model and it will ensure that they are reweighed with a constraint in mind. For this usecase, it would optimise;

\begin{equation} \label{eq1}
\begin{split}
\text{min} & \sum_i (y_i - (w_1 y_{1,i} + w_2 y_{2,i} + w_3 y_{3,i}))^2 \\
\text{subject to} & \\
& \text{ } w_1 + w_2 + w_3 = 1 \\ 
& \text{ } w_1 \geq 0, w_2 \geq 0, w_3 \geq 0
\end{split}
\end{equation}

The final positivity constraint is optional in our model. 

Here's an example usage.


```python
from sklearn.datasets import make_regression
import pandas as pd

X, y = make_regression(n_samples=1000, n_features=10, random_state=42)
df = pd.DataFrame(X)
```

We've turned the array into a dataframe so that we can apply the `ColumnSelector`.


```python
from sklearn.pipeline import Pipeline, FeatureUnion
from sklearn.linear_model import LinearRegression
from sklearn.decomposition import PCA
from sklearn.model_selection import GridSearchCV

from sklego.meta import EstimatorTransformer
from sklego.linear_model import ProbWeightRegression
from sklego.preprocessing import ColumnSelector

pipe = Pipeline([
    ('models', FeatureUnion([
        ('path1', Pipeline([
            ('select1', ColumnSelector([0,1,2,3,4])),
            ('pca', PCA(n_components=3)),
            ('linear', EstimatorTransformer(LinearRegression()))
        ])),
        ('path2', Pipeline([
            ('select2', ColumnSelector([5,6,7,8,9])),
            ('pca', PCA(n_components=2)),
            ('linear', EstimatorTransformer(LinearRegression()))
        ]))
    ])),
    ('prob_weight', ProbWeightRegression())
])

grid = GridSearchCV(estimator=pipe, param_grid={}, cv=3).fit(df, y)
```


```python
from sklearn import set_config
set_config(display='diagram')
grid
```




<style>div.sk-top-container {color: black;background-color: white;}div.sk-toggleable {background-color: white;}label.sk-toggleable__label {cursor: pointer;display: block;width: 100%;margin-bottom: 0;padding: 0.2em 0.3em;box-sizing: border-box;text-align: center;}div.sk-toggleable__content {max-height: 0;max-width: 0;overflow: hidden;text-align: left;background-color: #f0f8ff;}div.sk-toggleable__content pre {margin: 0.2em;color: black;border-radius: 0.25em;background-color: #f0f8ff;}input.sk-toggleable__control:checked~div.sk-toggleable__content {max-height: 200px;max-width: 100%;overflow: auto;}div.sk-estimator input.sk-toggleable__control:checked~label.sk-toggleable__label {background-color: #d4ebff;}div.sk-label input.sk-toggleable__control:checked~label.sk-toggleable__label {background-color: #d4ebff;}input.sk-hidden--visually {border: 0;clip: rect(1px 1px 1px 1px);clip: rect(1px, 1px, 1px, 1px);height: 1px;margin: -1px;overflow: hidden;padding: 0;position: absolute;width: 1px;}div.sk-estimator {font-family: monospace;background-color: #f0f8ff;margin: 0.25em 0.25em;border: 1px dotted black;border-radius: 0.25em;box-sizing: border-box;}div.sk-estimator:hover {background-color: #d4ebff;}div.sk-parallel-item::after {content: "";width: 100%;border-bottom: 1px solid gray;flex-grow: 1;}div.sk-label:hover label.sk-toggleable__label {background-color: #d4ebff;}div.sk-serial::before {content: "";position: absolute;border-left: 1px solid gray;box-sizing: border-box;top: 2em;bottom: 0;left: 50%;}div.sk-serial {display: flex;flex-direction: column;align-items: center;background-color: white;}div.sk-item {z-index: 1;}div.sk-parallel {display: flex;align-items: stretch;justify-content: center;background-color: white;}div.sk-parallel-item {display: flex;flex-direction: column;position: relative;background-color: white;}div.sk-parallel-item:first-child::after {align-self: flex-end;width: 50%;}div.sk-parallel-item:last-child::after {align-self: flex-start;width: 50%;}div.sk-parallel-item:only-child::after {width: 0;}div.sk-dashed-wrapped {border: 1px dashed gray;margin: 0.2em;box-sizing: border-box;padding-bottom: 0.1em;background-color: white;position: relative;}div.sk-label label {font-family: monospace;font-weight: bold;background-color: white;display: inline-block;line-height: 1.2em;}div.sk-label-container {position: relative;z-index: 2;text-align: center;}div.sk-container {display: inline-block;position: relative;}</style><div class="sk-top-container"><div class="sk-container"><div class="sk-item sk-dashed-wrapped"><div class="sk-label-container"><div class="sk-label sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="1f1cc5db-07d1-470d-aa5a-253d7d316602" type="checkbox" ><label class="sk-toggleable__label" for="1f1cc5db-07d1-470d-aa5a-253d7d316602">GridSearchCV</label><div class="sk-toggleable__content"><pre>GridSearchCV(cv=3,
             estimator=Pipeline(steps=[('models',
                                        FeatureUnion(transformer_list=[('path1',
                                                                        Pipeline(steps=[('select1',
                                                                                         ColumnSelector(columns=[0,
                                                                                                                 1,
                                                                                                                 2,
                                                                                                                 3,
                                                                                                                 4])),
                                                                                        ('pca',
                                                                                         PCA(n_components=3)),
                                                                                        ('linear',
                                                                                         EstimatorTransformer(estimator=LinearRegression()))])),
                                                                       ('path2',
                                                                        Pipeline(steps=[('select2',
                                                                                         ColumnSelector(columns=[5,
                                                                                                                 6,
                                                                                                                 7,
                                                                                                                 8,
                                                                                                                 9])),
                                                                                        ('pca',
                                                                                         PCA(n_components=2)),
                                                                                        ('linear',
                                                                                         EstimatorTransformer(estimator=LinearRegression()))]))])),
                                       ('prob_weight',
                                        ProbWeightRegression())]),
             param_grid={})</pre></div></div></div><div class="sk-parallel"><div class="sk-parallel-item"><div class="sk-item"><div class="sk-serial"><div class="sk-item"><div class="sk-serial"><div class="sk-item sk-dashed-wrapped"><div class="sk-label-container"><div class="sk-label sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="722adf51-b4eb-4aea-acce-ff4fd0835955" type="checkbox" ><label class="sk-toggleable__label" for="722adf51-b4eb-4aea-acce-ff4fd0835955">models: FeatureUnion</label><div class="sk-toggleable__content"><pre>FeatureUnion(transformer_list=[('path1',
                                Pipeline(steps=[('select1',
                                                 ColumnSelector(columns=[0, 1,
                                                                         2, 3,
                                                                         4])),
                                                ('pca', PCA(n_components=3)),
                                                ('linear',
                                                 EstimatorTransformer(estimator=LinearRegression()))])),
                               ('path2',
                                Pipeline(steps=[('select2',
                                                 ColumnSelector(columns=[5, 6,
                                                                         7, 8,
                                                                         9])),
                                                ('pca', PCA(n_components=2)),
                                                ('linear',
                                                 EstimatorTransformer(estimator=LinearRegression()))]))])</pre></div></div></div><div class="sk-parallel"><div class="sk-parallel-item"><div class="sk-item"><div class="sk-label-container"><div class="sk-label sk-toggleable"><label>path1</label></div></div><div class="sk-serial"><div class="sk-item"><div class="sk-serial"><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="ee9ceae5-59b7-4c88-9fbf-9505349fa887" type="checkbox" ><label class="sk-toggleable__label" for="ee9ceae5-59b7-4c88-9fbf-9505349fa887">ColumnSelector</label><div class="sk-toggleable__content"><pre>ColumnSelector(columns=[0, 1, 2, 3, 4])</pre></div></div></div><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="e0129048-6eb8-47bc-a8c0-ae9d96d63cd1" type="checkbox" ><label class="sk-toggleable__label" for="e0129048-6eb8-47bc-a8c0-ae9d96d63cd1">PCA</label><div class="sk-toggleable__content"><pre>PCA(n_components=3)</pre></div></div></div><div class="sk-item sk-dashed-wrapped"><div class="sk-label-container"><div class="sk-label sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="fe8293b5-d914-4a59-8c39-e2b4e18851c8" type="checkbox" ><label class="sk-toggleable__label" for="fe8293b5-d914-4a59-8c39-e2b4e18851c8">linear: EstimatorTransformer</label><div class="sk-toggleable__content"><pre>EstimatorTransformer(estimator=LinearRegression())</pre></div></div></div><div class="sk-parallel"><div class="sk-parallel-item"><div class="sk-item"><div class="sk-serial"><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="41c53c2a-d2c2-4241-b549-8f7653a0029f" type="checkbox" ><label class="sk-toggleable__label" for="41c53c2a-d2c2-4241-b549-8f7653a0029f">LinearRegression</label><div class="sk-toggleable__content"><pre>LinearRegression()</pre></div></div></div></div></div></div></div></div></div></div></div></div></div><div class="sk-parallel-item"><div class="sk-item"><div class="sk-label-container"><div class="sk-label sk-toggleable"><label>path2</label></div></div><div class="sk-serial"><div class="sk-item"><div class="sk-serial"><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="86b9c93e-12be-4fa6-865f-bb874ff99129" type="checkbox" ><label class="sk-toggleable__label" for="86b9c93e-12be-4fa6-865f-bb874ff99129">ColumnSelector</label><div class="sk-toggleable__content"><pre>ColumnSelector(columns=[5, 6, 7, 8, 9])</pre></div></div></div><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="ec46eb32-4807-466f-a47b-dba193ff5315" type="checkbox" ><label class="sk-toggleable__label" for="ec46eb32-4807-466f-a47b-dba193ff5315">PCA</label><div class="sk-toggleable__content"><pre>PCA(n_components=2)</pre></div></div></div><div class="sk-item sk-dashed-wrapped"><div class="sk-label-container"><div class="sk-label sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="5efe2360-3f3b-4c67-bc8f-e2619f06d229" type="checkbox" ><label class="sk-toggleable__label" for="5efe2360-3f3b-4c67-bc8f-e2619f06d229">linear: EstimatorTransformer</label><div class="sk-toggleable__content"><pre>EstimatorTransformer(estimator=LinearRegression())</pre></div></div></div><div class="sk-parallel"><div class="sk-parallel-item"><div class="sk-item"><div class="sk-serial"><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="6ceb127c-c870-4944-b6c3-5f47c0b5af1b" type="checkbox" ><label class="sk-toggleable__label" for="6ceb127c-c870-4944-b6c3-5f47c0b5af1b">LinearRegression</label><div class="sk-toggleable__content"><pre>LinearRegression()</pre></div></div></div></div></div></div></div></div></div></div></div></div></div></div></div><div class="sk-item"><div class="sk-estimator sk-toggleable"><input class="sk-toggleable__control sk-hidden--visually" id="8d1aaabb-dda3-4fe0-99dc-c4b60597b332" type="checkbox" ><label class="sk-toggleable__label" for="8d1aaabb-dda3-4fe0-99dc-c4b60597b332">ProbWeightRegression</label><div class="sk-toggleable__content"><pre>ProbWeightRegression()</pre></div></div></div></div></div></div></div></div></div></div></div></div>



You can see that the `ProbWeightRegression` indeeds sums to one. 


```python
grid.best_estimator_[1].coefs_
```




    array([0.03102466, 0.96897535])



## Least Absolute Deviation Regression

Imagine that you have a dataset with some outliers.


```python
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(0)
X = np.linspace(0, 1, 20)
y = 3*X + 1 + 0.5*np.random.randn(20)
X = X.reshape(-1, 1)

y[10] = 8
y[15] = 15

plt.figure(figsize=(16, 4))
plt.scatter(X, y)
```

A simple linear regression will not do a good job since it is distracted by the outliers. That is because it optimizes the mean squared error
$$ \sum_i \left(y_i-\textrm{model}(x_i)\right)^2 $$

which penalizes a few large errors more than many tiny errors. For example, if y-model(x) = 4 for some single observation, the MSE here is 16. If there are two observations with y_1 - model(x_1) = 2 and y_2 - model(x_2) = 2, the MSE is 8 in total, which is less than for one larger error. Note that the sum of the errors is the same in both cases.

Hence, linear regression does the following:


```python
from sklearn.linear_model import LinearRegression

x = np.array([0, 1]).reshape(-1, 1)
plt.figure(figsize=(16, 4))
plt.scatter(X, y)
plt.plot(x, LinearRegression().fit(X, y).predict(x), 'r')
```

By changing the loss function to the mean absolute deviation 
$$ \sum_i \left|y_i-\textrm{model}(x_i)\right|\enspace, $$

we can let the model put the same focus on each error. This yields the least absolute deviation (LAD) regression that tries to agree with the majority of the points.


```python
from sklearn.linear_model import LinearRegression
from sklego.linear_model import LADRegression

x = np.array([0, 1]).reshape(-1, 1)
plt.figure(figsize=(16, 4))
plt.scatter(X, y)
plt.plot(x, LinearRegression().fit(X, y).predict(x), 'r')
plt.plot(x, LADRegression().fit(X, y).predict(x), 'g')
```

### See also
scikit-learn tackles this problem by offering a variety of [robust regressors](https://scikit-learn.org/stable/auto_examples/linear_model/plot_robust_fit.html). Many of them use an indirect approach to reduce the effect of outliers. [RANSAC](https://en.wikipedia.org/wiki/Random_sample_consensus), for example, samples random points from the dataset until it consists of all inliers.
The closest thing to LADRegression that scikit-learn offers is the [HuberRegressor](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.HuberRegressor.html) with a loss function that is partly a squared and partly an absolute error. However, it is more complicated and requires hyperparameter tuning to unleash its full potential.


## QuantileRegression

This is an extension of the LADRegression (see above). While the LADRegression outputs a line that over- and underestimates around 50% of the data, the Quantileregression yields lines with different over- and underestimation shares. This can be used for creating simple confidence intervals around predictions. As an example, consider the following:
1. Create a QuantileRegression with quantile=0.1,
2. create a QuantileRegression with quantile=0.9,

and around 80% of the data is between these two lines.


```python
from skbonus.linear_model import QuantileRegression
import numpy as np
import matplotlib.pyplot as plt

np.random.seed(123)
X = np.arange(100).reshape(-1, 1)
y = 2*X.ravel() + X.ravel()*np.random.standard_cauchy(100)

q_10 = QuantileRegression(quantile=0.1).fit(X, y)
q_90 = QuantileRegression(quantile=0.9).fit(X, y)
lad = QuantileRegression().fit(X, y)

plt.plot(X, y)
plt.plot(X, lad.predict(X))
plt.fill_between(X.ravel(), q_10.predict(X), q_90.predict(X), alpha=0.33, color='orange')
```


```python

```

98fd62369c950f41214b6d4f0532a642b53adda3
ipynb
Jupyter Notebook
doc/linear-models.ipynb
lahdjirayhan/scikit-lego
5dd145df796c4d254cd505727c9db01484ebc39c
MIT
2019-03-01T21:35:53.000Z
2022-03-30T11:22:46.000Z
doc/linear-models.ipynb
lahdjirayhan/scikit-lego
5dd145df796c4d254cd505727c9db01484ebc39c
MIT
2019-02-27T10:38:53.000Z
2022-03-31T07:22:24.000Z
doc/linear-models.ipynb
lahdjirayhan/scikit-lego
5dd145df796c4d254cd505727c9db01484ebc39c
MIT
2019-03-01T19:34:37.000Z
2022-03-30T14:10:29.000Z
Qwen/Qwen-72B
1. YES
2. YES