ใช้ Bayesian ทำ Linear Regression แล้วหน้าตาออกมาเป็นยังไง

4 min readAug 19, 2020

ได้ยินมานานแล้วสำหรับ Bayesian probability แต่ไม่เคยได้ศึกษาจริงจังเลย เพราะที่ผ่านมาส่วนมากพยายามยึดโยงสิ่งต่าง ๆ เข้ากับงานที่ทำอยู่ แต่ยังหาจุดร่วมระหว่างงานที่ทำกับ Bayesian ไม่ได้สักที แต่วันนี้มาขอลองทำแล้วก็หยิบงานที่เคยทำ เช่น Linear regression model (เข้าใจง่าย) เอามาเป็นเคสตัวอย่าง

Frequentist vs Bayesian

เราควรมาเริ่มทำความรู้จักก่อนว่า Bayesian คือใครหรือคืออะไร ดังนั้นเราควรรู้ก่อนว่า จริง ๆ แล้ว Linear regression model ที่รู้จักมาตั้งแต่มัธยม เป็นการทำโมเดลแบบ Frequentist’s view โดยผลลัพธ์จากโมเดล เกิดจากการหาค่าที่ทำให้ Minimize the residual sum of squares (RSS) ดังนั้นค่า “คงที่” หรือที่เรียกว่า Deterministic คือมีค่าเป็นจุดหนึ่งจุดบนเส้นตรง

ดังนั้นเมื่อแทนค่า Features ลงไปในสมการที่มีค่าคงที่ (Co-efficient) โมเดลก็ทำการ Estimate ค่าออกมา ตามสมการข้างล่าง โดยที่ตัวแปรเดียวหรือหลายตัวแปร ก็มีหลักการทำงานเหมือนกัน

ในขณะที่ทฤษฏี Bayesian ไม่ได้ทำการหาค่าแบบ Deterministic (คงที่) ลักษณะการคำนวณคือการประมาณค่าจาก Probabilistic Distribution หรือค่าของ Weight เกิดจาก Gaussian Distribution (Normal Distribution)

The response, y, is not estimated as a single value, but is assumed to be drawn from a probability distribution.

เนื่องจากทฤษฏีเป็นแบบนี้ ทำให้การสร้างโมเดลเลยเปลี่ยนไป เพราะการหาค่า Y ไม่ได้เกิดจากการหาค่าคงที่เพียง 1 ค่าเหมือนเดิมเหมือนแล้ว แต่ค่าของ Y ขึ้นกับความเป็นไปได้ทั้งหมดของ X และ Weight

Install

เนื่องจากวันนี้เราต้องการทำ Linear model โดยใช้ Bayesian จึงจำเป็นต้องติดตั้ง Library เพิ่มเติมที่มีชื่อว่า PyMC3 สามารถอ่าน Documentation เพิ่มเติมได้จากลิงค์ข้างล่าง

PyMC3 Documentation - PyMC3 3.8 documentation

Sometimes an unknown parameter or variable in a model is not a scalar value or a fixed-length vector, but a function. A…

docs.pymc.io

เราต้องการ Visualize ผลลัพธ์ของ Bayesian model ดังนั้นจำเป็นต้องติดตั้ง Library เพิ่มเติมคือ ArviZ เพื่อทำการ Plot กราฟต่าง ๆ ออกมาดูเพื่อให้เห็นภาพที่ชัดเจนกว่าเดิม

ArviZ: Exploratory analysis of Bayesian models - ArviZ dev documentation

ArviZ is a Python package for exploratory analysis of Bayesian models. Includes functions for posterior analysis, data…

arviz-devs.github.io

เนื่องจาก Library อาจมี Bug ในเวอร์ชั่นใหม่ จึงขอแนะนำให้ติดตั้งโดยระบุเวอร์ชั่นที่ต้องการใช้งาน โดยวันนี้ขอใช้งานเวอร์ชั่นตามที่เห็นข้างล่าง

!pip install arviz==0.6.1
!pip install pymc3==3.8

Code

ความยึดโยงเดียวที่นึกออกคือ การนำโมเดล Linear regression ที่เคยทำไว้แล้ว มาเปรียบเทียบกับการทำโมเดลด้วยทฤษฎี Bayesian ว่าผลที่ได้จะมีความแตกต่างกันมากน้อยแค่ไหน (หรืออาจไม่ต่างเลย) โดยโมเดลที่รันด้วย Linear regression ปกติได้ผลลัพธ์ตามด้านล่าง

R_Square: 0.8616087471847003
MSE: 0.018211503105028705
                         OLS Regression Results                             ==================================================================== Dep. Variable:            y   R-squared:                       0.862 Model:                  OLS   Adj. R-squared:                  0.853 Method:       Least Squares   F-statistic:                     97.54 Date:      Mon, 10 Aug 2020   Prob (F-statistic):           3.39e-20 Time:              13:42:48   Log-Likelihood:                 31.862 No. Observations:        51   AIC:                            -55.72 Df Residuals:            47   BIC:                            -48.00 Df Model:                 3                                          Covariance Type:  nonrobust                                          ====================================================================
          coef   std err        t       P>|t|      [0.025     0.975] --------------------------------------------------------------------const  -3.6165     0.022 -161.740       0.000      -3.662     -3.572
x1     23.2338     9.115    2.549       0.014       4.897     41.570
x2     -0.1377     0.048   -2.849       0.006      -0.235     -0.040
x3      3.1357    0.322     9.728       0.000       2.487      3.784 ====================================================================
Omnibus:              7.794   Durbin-Watson:                   0.773 Prob(Omnibus):        0.020   Jarque-Bera (JB):                7.187 Skew:                -0.697   Prob(JB):                       0.0275 Kurtosis:             4.200   Cond. No.                         483. ====================================================================

มาลองใช้ Bayesian สร้างโมเดลโดยใช้ตัวแปรเดียวกัน

การใช้งาน Library PyMC3 ในการ Simulate Bayesian model สามารถทำได้ด้วยการ Define model ด้วยคำสั่ง with… as … : แล้วสร้างตัวแปรมารับ 1 ตัว ซึ่งวันนี้ขอเรียกว่า normal_model จากนั้นต้องกำหนด Distribution ของโมเดลว่าต้องการให้อยู่ใน Curve แบบไหน ซึ่งวันนี้กำหนดเป็น Normal distribution และเอาตัวแปร family มารับค่าไว้

ต่อมาเป็นการกำหนด Formula ให้กับตัวโมเดลว่า Effect ของ x เป็นอย่างไรกับ y ซึ่งสามารถเขียนเหมือน Formula ของ Linear model ได้ว่า y ~ x1 + x2 + x3 จากนั้นก็ใส่ Input (ในที่นี้เป็น DataFrame) ที่มีข้อมูลทั้งตัวแปร x และ y เข้าไป และประกาศfamily ให้เป็น Normal distribution ตามที่เราสร้างไว้

ใช้ pm.sample() ในการสร้าง (Simulation) โมเดลออกมา โดยให้ Markov Chain Monte Carlo (MCMC) เป็นตัวช่วย Simulate (Simulator) ถ้าทำทุกอย่างถูกต้อง จะมีหน้าต่าง Processing ขึ้นมาเหมือนข้างล่าง

Auto-assigning NUTS sampler...
Initializing NUTS using jitter+adapt_diag...
Multiprocess sampling (2 chains in 4 jobs)
NUTS: [sd, x3, x2, x1, Intercept]
Sampling 2 chains, 0 divergences: 100%|██████████| 5000/5000
[00:09<00:00, 550.13draws/s]

Result

เมื่อรันคำสั่ง pm.sample() เสร็จแล้ว Bayesian ก็ถูกสร้างเอาไว้เป็นที่เรียบร้อย ซึ่งสามารถทำการ Plot เพื่อดูผลลัพธ์ที่เกิดขึ้นได้

เนื่องจากผลลัพธ์ออกมาในรูปแบบ Normal distribution ดังนั้นเราต้องหาค่ากลาง (Mean) ออกมา เพื่อใช้ในการ Estimate ผลลัพธ์

ใช้คำสั่ง pm.summary() ได้เป็น Dataframe ที่มีการสรุปค่าจากโมเดลทั้งหมด เมื่อลองเทียบกับผลลัพธ์จาก Linear regression ปกติที่เคยทำ พบว่าค่าที่ได้ออกมามีความใกล้เคียงกันมาก

Prediction

ทดลองใช้งาน Bayesian linear model ตัวอย่างวันนี้คือทดลอง Simulate ทั้งหมด 100 Samples

100%|██████████| 100/100 [00:00<00:00, 115.23it/s]

ตัวแปร ppc (ย่อมาจาก Posterior Predictive Check) ออกมาเป็น Dictionary ที่เก็บค่าทั้งหมด 100 ตัวอย่างตามจำนวน Observation ที่เรา Train โมเดลเอาไว้ เพื่อให้เห็นภาพชัดเจนขึ้น ขอ Plot เส้นทั้งหมด 100 เส้นออกมา

กราฟข้างบนคือความเป็นไปได้ทั้งหมด ตาม Normal distribution ที่เราสุ่มค่าเข้าไปในโมเดล

เส้นสีแดงคือค่าเฉลี่ยของ Predicted y ที่เกิดจากการ Simulate ทั้งหมด

จุด (Area) ที่เป็นสีแดงคือค่าที่ออกมาจากโมเดลทั้งหมด 100 ชุด และเส้นสีน้ำเงินคือค่าจริง (y) ที่เรากำลังพยายาม Predict อยู่

ลองมาเทียบกับผลลัพธ์จาก Linear regression ปกติดูบ้าง

จากกราฟข้างบนก็เห็นแล้วว่า Pattern ของโมเดลยังเป็นไปในทิศทางเดียวกันอยู่ ซึ่งก็ Support ได้จากค่า Mean ของ Estimate ต่าง ๆ ขอจบการลองทำ Bayesian ไว้แต่เพียงเท่านี้