رگرسیون چندجمله‌ای در پایتون

در این نوشتار، مروری بر الگوریتم رگرسیون چندجمله‌ای خواهیم داشت. این الگوریتم با تغییر تابع فرضیه و افزودن چند ویژگی جدید به ورودی‌ها، برای برازش داده‌های غیرخطی به کار می‌رود. رگرسیون چندجمله‌ای نسخه‌ای از رگرسیون خطی استاندارد است.

در بخش اول که مربوط به معرفی نمادهاست، نسبت به قسمت‌های قبلی، یک متغیر جدید به نام degrees را معرفی می‌کنیم. در بخش بعدی، تعریفی از الگوریتم ارائه می‌دهیم. در انتها، با استفاده از NumPy و Matplotlib پایتون، پیش‌بینی‌های انجام‌شده را مصورسازی و نمره‌ r2 را محاسبه خواهیم کرد.

نمادها

• n: تعداد ویژگی‌ها
• m: تعداد نمونه‌های آموزشی
• X: ماتریس ورودی با اندازه‌ mxn
• y: بردار مقدار حقیقی/هدف با اندازه‌ m
• x(i), y(i): iمین نمونه‌ آموزشی؛ x(i)، nبعدی و y(i) یک عدد حقیقی است.
• degrees: لیست؛ ویژگی X^(value) به ورودی افزوده می‌شود، value از داخل لیست انتخاب می‌شود (توضیحات بیشتر در قسمت‌های بعدی).
• w: وزن‌ها (پارامترها) با شکل nx1
• b: سوگیری (پارامتر)، یک عدد حقیقی که می‌تواند انتشار همگانی یابد.
• y_hat: فرضیه (ضرب ماتریسی w و x به‌علاوه‌ b) یا X + b
• تابع زیان: تابع زیان MSE یا خطای میانگین مجذورات (y_hat-y)²؛ برای پیدا کردن پارامترهای w و b باید این تابع را به حداقل رساند.

رگرسیون چندجمله‌ای

برای درک بهتر رگرسیون چندجمله‌ای، از دیتاست پایین استفاده می‌کنیم. محور x  نشان‌دهنده‌ داده‌های ورودی x و محور y نشان‌دهنده‌ مقادیر حقیقی/هدف y است. در این دیتاست 1000 نمونه‌ (m) و 1 ویژگی (n) وجود دارد.

import numpy as np
import matplotlib.pyplot as pltnp.random.seed(42)
X = np.random.rand(1000,1)
y = 5*((X)**(2)) + np.random.rand(1000,1)

با استفاده از رگرسیون خطی، تنها می‌توان یک خط صاف روی این داده‌ها برازش داد (خط آبی در تصویر پایین) و فرضیه چنین شکلی خواهد داشت: w1.X + b (w1 به جای w).

اما همان‌طور که در تصویر مشاهده می‌کنید، از آنجایی که داده‌ها خطی نیستند، خط قرمزرنگ بهتر روی داده‌ها برازش یافته است.

سؤالی که هنگام برازش مدل روی این داده‌ها مطرح می‌شود این است که «از کدام ویژگی‌ها باید استفاده کرد؟». آیا باید از یک خط صاف برای برازش داده‌ها استفاده کرد یا از آنجایی که داده‌ها به شکل یک تابع درجه‌ دو هستند، باید از فرضیه‌ای به شکل b + w1.X + w2.X²  استفاده کرد؟

شاید هم به خاطر شکل داده‌ها، بهتر است (تابع) فرضیه یک تابع جذری باشد: b + w1.X + w2.(X)^0.5. در واقع، داده‌ها می‌توانند هر درجه‌ای باشند، شما می‌توانید ویژگی‌های موردنظر را انتخاب کرده و تابع فرضیه را مطابق با آن تغییر دهید.

برای پیاده‌سازی این تابع، تنها کاری که باید انجام دهیم، تعریف اولین ویژگی x1 به‌عنوان x و دومین ویژگی x2 به‌عنوان X² است. راه دیگر این است که x1 را به‌صورت X^0.5 تعریف کنید؛ اینکه کدام راه را انتخاب می‌کنید، بستگی به ویژگی‌هایی دارد که قصد دارید به کار ببرید. با تعریف ویژگی جدید x2 به‌صورت X² یا X^0.5 می‌توان از همان مکانیزمی که در رگرسیون خطی مشاهده کردیم، برای برازش این داده‌های غیرخطی نیز استفاده کنیم.

نکته‌ مهم این است که فرضیه‌ ما همچنان خطی است، چون X² یا X^0.5  صرفاً ویژگی هستند. بااین‌حال، با استفاده از آن می‌توانیم برازش غیرخطی روی داده‌ها انجام دهیم.

تنها کاری که باید انجام دهیم، تغییر ورودی x، یعنی افزودن ویژگی‌های درجه‌ اول، دوم، و… به آن است؛ ما در این مثال، یک ویژگی X² (درجه‌ دوم) را به فرضیه اضافه کرده‌ایم.

اگر ورودی را از طریق افزودن ویژگی‌های جدید تغییر دهید، فرضیه خودبه‌خود تغییر می‌کند، زیرا h(x) = w.X +b وw  برداری با اندازه‌ n (تعداد ویژگی‌ها) است.

شکل پایین ده نمونه‌ اول دیتاست را بعد از افزودن ویژگی جدید X² به ورودی، نشان می‌دهد:

هر تعداد ویژگی که بخواهیم می‌توانیم به داده‌ها اضافه کنیم، به شرطی که توانی از ویژگی موجود باشد.

الگوریتم

• تغییر تابع فرضیه رگرسیون خطی متناسب با داده‌های موجود
• افزودن ویژگی‌های جدید با توان (درجه‌) بالاتر به ورودی‌ها
• اجرای الگوریتم گرادیان کاهشی/mini-batch روی ورودی‌های تغییریافته به‌منظور پیدا کردن پارامترهای b و w (سوگیری و وزن)

اجرای رگرسیون چندجمله‌ای

تابع زیان

ابتدا تابع زیان یعنی MSE را تعریف می‌کنیم:

(y_hat-y)²  که در آن y_hat فرضیه w.X + b است.

def loss(y, y_hat):
    
    # y --> true/target value.
    # y_hat --> hypothesis
    
    #Calculating loss.
    loss = np.mean((y_hat - y)**2)
    return loss

تابعی برای محاسبه‌ گرادیان‌ها

حال برای محاسبه‌ مشتق‌های جزئی (گرادیان‌) تابع زیان نسبت به متغیرهای w و b، یک تابع می‌نویسیم؛ به کامنت‌ها دقت کنید.

# Calulating gradient of loss w.r.t parameters(weights and bias).def gradients(X, y, y_hat):
    
    # X --> Input.
    # y --> true/target value.
    # y_hat --> hypothesis
    # w --> weights (parameter).
    # b --> bias (parameter).
    
    # m-> number of training examples.
    m = X.shape[0]
    
    # Gradient of loss w.r.t weights.
    dw = (1/m)*np.dot(X.T, (y_hat - y))
    
    # Gradient of loss w.r.t bias.
    db = (1/m)*np.sum((y_hat - y)) 
    
    return dw, db

تابعی برای افزودن ویژگی به داده‌های ورودی

در این قسمت تابعی می‌نویسیم که اجازه می‌دهد هر تعداد ویژگی که می‌خواهیم را به ورودی اضافه کنیم. بدین منظور یک متغیر به نام degrees تعریف می‌کنیم؛ این متغیر یک لیست پایتونی است.

هر مقداری که به این لیست وارد کنیم، یک ویژگی X^(value) به ورودی اضافه می‌شود. برای مثال، اگر 2 و 3 را وارد لیست کنیم، ویژگی‌های X² و X³ را به ورودی اضافه کرده‌ایم.

def x_transform(X, degrees):
    
    # X --> Input.
    # degrees --> A list, We add X^(value) feature to the input
    #             where value is one of the values in the list.
    
    # making a copy of X.
    t = X.copy()
    
    # Appending columns of higher degrees to X.
    for i in degrees:
        X = np.append(X, t**i, axis=1)
            
    return X

تابع آموزش

تابع آموزش، وزن‌ها، سوگیری‌ها و حلقه‌ آموزشی را با گرادیان کاهشی mini-batch تعریف می‌کند.

def train(X, y, bs, degrees, epochs, lr):
    
    # X --> Input.
    # y --> true/target value.
    # bs --> Batch Size.
    # epochs --> Number of iterations.
    # degrees --> A list, We add X^(value) feature to the input
    #             where value is one of the values in the list.
    # lr --> Learning rate.
    
    # Adding features to input X.
    x = x_transform(X, degrees)
    
    # m-> number of training examples
    # n-> number of features 
    m, n = x.shape
    
    # Initializing weights and bias to zeros.
    w = np.zeros((n,1))
    b = 0
    
    # Reshaping y.
    y = y.reshape(m,1)
    
    # Empty list to store losses.
    losses = []
    
    # Training loop.
    for epoch in range(epochs):
        for i in range((m-1)//bs + 1):
            
            # Defining batches.
            start_i = i*bs
            end_i = start_i + bs
            xb = x[start_i:end_i]
            yb = y[start_i:end_i]
            
            # Calculating hypothesis
            y_hat = np.dot(xb, w) + b
            
            # Getting the gradients of loss w.r.t parameters.
            dw, db = gradients(xb, yb, y_hat)
            
            # Updating the parameters.
            w -= lr*dw
            b -= lr*db
        
        # Calculating loss and appending it in the list.
        l = loss(y, np.dot(x, w) + b)
        losses.append(l)
        
    # returning weights, bias and losses(List).
    return w, b, losses

تابع پیش‌بینی

# Predicting function.def predict(X, w, b, degrees):
    
    # X --> Input.
    # w --> weights (parameter).
    # b --> bias (parameter).
    degrees --> A list, We add X^(value) feature to the input
    #             where value is one of the values in the list.
    
    # Adding degrees to input X.
    x1 = x_transform(X, degrees)
    
    # Returning predictions.
    return np.dot(x1, w) + b

آموزش و مصورسازی پیش‌بینی‌ها

آموزش داده‌ها با استفاده از تابع train انجام می‌گیرد.

اینجا فقط 2 را وارد لیست degrees کرده‌ایم. شما هر عدد یا اعدادی را که بخواهید، می‌توانید وارد کنید.

w, b, l = train(X, y, bs=100, degrees=[2], epochs=1000,
                lr=0.01)# Plotting
fig = plt.figure(figsize=(8,6))
plt.plot(X, y, 'y.')
plt.plot(X, predict(X, w, b, [2]), 'r.')
plt.legend(["Data", "Polynomial predictions"])
plt.xlabel('X - Input')
plt.ylabel('y - target / true')
plt.title('Polynomial Regression')
plt.show()

از آنجایی که تابع زیان مربوط به هر دور آموزش را در لیست losses جمع‌آوری کرده‌ایم، اکنون می‌توانیم نمودار این لیست را به ازای دورها (epochs) رسم کرده و ببینیم آیا تابع زیان در طی آموزش کاهش یافته است یا خیر.

نمره r2

محاسبه‌ نمره‌ r2 برای پیش‌بینی‌ها به ارزیابی عملکرد کمک می‌کند.

def r2_score(y, y_hat):
    return 1 - (np.sum((np.array(y_hat)-np.array(y))**2)/
                np.sum((np.array(y)-np.mean(np.array(y)))**2))r2_score(y_train, predict(x_train, w, b, [2]))
>>0.9541881152879292

نمره r2 که با استفاده از این کد به دست آمده، حاکی از عملکرد خوب مدل است.