[np-ml] Linear Regression

合集 - numpy-ml(2)

1.[np-ml] Ridge Regression2024-09-14

2.[np-ml] Linear Regression2024-09-14

收起

算法描述

A weighted linear least-squares regression model.

Notes
In weighted linear least-squares regression, a real-valued target vector, $\mathbf{y}$, is modeled as a linear combination of covariates $\mathbf{X}$, and model coefficients, $\beta$:

$$y_i={\beta }^{\mathrm{\top }}{\mathbf{x}}_i+{ϵ}_i$$

where ${ϵ}_i\sim \mathcal{N}(0,{\sigma }_i^2)$ is the error term associated with $i$-th example, and ${\sigma }_i^2$ is the variance of the corresponding example.

Under the model, the maximum-likelihood estimate for the regression coefficients, $\beta$, is:

$$\hat{\beta }={\mathrm{\Sigma }}^{-1}{\mathbf{X}}^{\mathrm{\top }}\mathbf{Wy}$$

where ${\mathrm{\Sigma }}^{-1}=({\mathbf{X}}^{\mathrm{\top }}\mathbf{WX}{)}^{-1}$ and $\mathbf{W}$ is a diagonal matrix of weights, with each entry inversely proportional to the variance of the corresponding measurement.
When $\mathbf{W}$ is the identity matrix the examples are weighted equally and the model reduces to standard linear least squares.

算法实现

直接按照解析解进行更新：

import numpy as np
 
 
def fit(X, y, fit_intercept=True, weights=None):
    N = X.shape[0]
 
    weights = np.ones(N) if weights is None else np.atleast_1d(weights)
    weights = np.squeeze(weights) if weights.size > 1 else weights
    err_str = f"`weights` must have shape ({N},), but got {weights.shape}"
    assert weights.shape == (N,), err_str
 
    # scale X and y by the weight associated with each examples
    W = np.diag(np.sqrt(weights))
    X, y = W @ X, W @ y
 
    # convert X to a design matrix if we are fitting an intercept
    if fit_intercept:
        X = np.c_[np.sqrt(weights), X]
 
    sigma_inv = np.linalg.pinv(X.T @ X)
    beta = np.atleast_2d(sigma_inv @ X.T @ y)
 
    return sigma_inv, beta

---

SGD方式进行模型参数更新：

import numpy as np
from numpy import ndarray
 
 
def main(X: ndarray, y: ndarray, num_iteration=1000, learning_rate=1e-2, weights=None):
    _, input_dim = X.shape
 
    # initialize weights
    weights = weights or np.random.randn(input_dim, 1)
 
    for iteration in range(num_iteration):
        # (batch_size, 1)
        y_pred = X @ weights
        # (batch_size, 1)
        diff = y_pred - y
        # (input_dim, 1)
        grad = X.T @ diff
 
        weights -= learning_rate * grad
 
    return weights
 
 
def test_main():
    A = np.array([[1, 1, 0]]).T
    X = np.random.rand(5, 3)
    y = X @ A
    weights = main(X, y, num_iteration=10000, learning_rate=1e-2)
    assert np.allclose(A, weights, rtol=1e-4), f"{weights=}\n but expect {A=}"
 
 
if __name__ == "__main__":
    test_main()

参考资料

code
doc