Forecast ราคาน้ำมัน จากข้อมูลย้อนหลังด้วย SARIMA Model

เคยเขียนเรื่อง ARIMA Model ไปทั้งหมด 6 ตอนแล้ว ลองมาเพิ่ม Assumption ให้โมเดลโดยยังคง Concept เดิมไว้ ซึ่ง Assumption ที่พูดถึงคือ Seasonal Auto Regressive Integrated Moving Average (SARIMA) หรือ Seasonal ARIMA

ก่อนเข้าสู่รายละเอียดของ SARIMA Model ขอย้อนกลับไปในเรื่องราวของการ ARIMA Model กันก่อน อย่างที่รู้กันดีอยู่แล้วว่า ARIMA Model คือการพยายามทำนายข้อมูลในอนาคต ด้วยข้อมูลของตัวเองจากอดีตที่เกิดขึ้นแล้ว ซึ่ง ARIMA Model พิจารณาแนวโน้ม หรือ Trend เพื่อที่จะหา Movement ที่อาจเกิดขึ้นในอนาคต แต่เนื่องจากข้อมูลอาจไม่ได้มีเพียงแค่ Trend ที่อยู่ในนั้น แต่อาจมี Seasonal effects ปะปนอยู่ด้วย จึงเป็นที่มาของ Seasonal ARIMA

SARIMA Model

รูปแบบของ SARIMA Model เขียนได้ว่า SARIMA(p,d,q)(P,D,Q)m ซึ่ง p, d, q เคยรู้จักกันมาก่อนหน้านี้แล้วจาก ARIMA Model ดังนั้นยังเหลืออีก 4 ตัวที่ต้องทำความรู้จักกันใหม่

A seasonal ARIMA model is formed by including additional seasonal terms in the ARIMA […] The seasonal part of the model consists of terms that are very similar to the non-seasonal components of the model, but they involve backshifts of the seasonal period.

ซึ่งจริง ๆ แล้ว P, D, Q มีความหมายเหมือนเดิมตามที่เคยรู้จักก่อนหน้านี้ ไม่ว่าจะเป็น P ที่เป็นตัวแทนของ Lag time, D ที่เป็นตัวแทนของ Integrated และ Q ที่เป็นตัวแทนของ Lag time of error แต่ทั้งหมดนี้เกิดขึ้นบน Seasonal แทนที่จะเป็น Trend ตัวที่ต้องสนใจต่อมาคือ “m” หรืออาจเขียนเป็น “s” ให้แทนว่า Seasonal ไปเลยก็ได้ สังเกตว่า m หรือ s อยู่นอกวงเล็บ ดังนั้นหมายความว่า Effect ของ Seasonal ต้องคูณเข้าไปทั้ง P, D และ Q ขอยกตัวอย่างให้เห็นภาพมากขึ้น

ตัวอย่างเช่น (P,D,Q)m = (1,1,0)12 หมายความว่า Effect ของ P, D, Q เกิดขึ้นจาก 12 periods ก่อนหน้า ดังนั้นถ้าข้อมูลในโมเดลเป็น Monthly basis สิ่งที่เกิด ณ ตอนนี้เป็นผลจาก 12 เดือนที่แล้ว ซึ่งคือ 1 ปีก่อนหน้านั่นเอง

Data and Assumption

โดยทั่วไปแล้ว Autoregressive model เหมาะกับการใช้ทำนายอนาคตกับสิ่งที่หาความสัมพันธ์ได้ยาก หรือความสัมพันธ์ที่มีผลอาจไม่ได้อยู่ในรูปแบบตัวเลข ซึ่งยากต่อการสร้างโมเดลคณิตศาสตร์ ครั้งนี้ขอใช้ข้อมูลราคาน้ำมัน Crude Oil Prices: Brent — Europe ในการทดลองสร้างโมเดลขึ้นมา ใช้ข้อมูลย้อนหลังทั้งหมด 5 ปี โดยเริ่มจาก 28 August 2015 ถึง 28 August 2020

Crude Oil Prices: Brent - Europe

กำหนด Assumption ให้ Seasonal effect เกิดขึ้น 1 อาทิตย์ที่ส่งผลถึงปัจจุบัน เนื่องจากข้อมูลเป็น Daily basis จึงกำหนดให้ Seasonal มีค่าเท่ากับ 7 ตรงนี้สามารถเปลี่ยน Effect เป็น 1 เดือน, 3 เดือน, 1 ไตรมาส หรือ 1 ปีก็ทำได้เช่นกัน (30, 90, 180, 365)

ตัวแปรอื่น ๆ เช่น p, d, q และ P, D, Q ขอกำหนดให้อยู่ช่วงใน 0–2 (0,1,2) เพราะโดยทั่วไปแล้ว Movement ของข้อมูลจะอยู่ในช่วงนี้

Code

เนื่องจากต้องการใช้ Grid-search ให้การหาค่าที่ Optimize ที่สุดจากค่า AIC ทำให้ Require run-time ที่ค่อนข้างนาน เนื่องจาก Combinations ที่เยอะขึ้น จึงขอใช้ Colab ในการทำโมเดล

Library ที่จำเป็นต่อการใช้งานคือ SARIMAX ซึ่งเป็นหนึ่งใน API ของ Statsmodels สามารถอ่านรายละเอียดเพิ่มเติมของฟังก์ชั่นนี้ได้ที่ Link ด้านล่าง เริ่มอ่าน Dataset เข้ามาในตัวแปรตามปกติของการเริ่มต้นทำโมเดล ขอ .set_index() ด้วย Column DATE ที่ parse_dates เพราะต้องการให้เหลือเฉพาะ Column ที่ต้องการใช้งานคือ price

statsmodels.tsa.statespace.sarimax.SARIMAX - statsmodels

ทำ Data assessment แล้วเจอว่าข้อมูลราคาน้ำมันที่อ่านเข้ามา ถูกเก็บค่าไว้เป็น String ดังนั้นต้องแปลงค่าก่อนถึงเริ่มการคำนวณได้ นอกจากนี้ยังมี Missing value ที่ถูกเก็บด้วย . อยู่ประมาณ 30 วัน ดังนั้นก่อนเริ่มทำโมเดลจึงจำเป็นต้องแก้ไขสิ่งผิดพลาดเหล่านี้ก่อน

เนื่องจาก String ที่มีค่า Missing value ที่เก็บเป็น . ไม่สามารถเปลี่ยนแปลงเป็นค่า Float ได้ทันที ดังนั้นต้องทำการแทนค่า . ด้วย numpy.nan โดยใช้ .replace() จากนั้นจึงสามารถเปลี่ยนทั้ง Column price เป็น numpy.float64 โดยใช้ .astype() สุดท้ายเป็นแทนค่า Missing value ด้วยตัวเลข ซึ่งขอให้ค่าเฉลี่ยเพื่อแทนค่าลงไป สามารถใช้ .fillna(.mean()) เมื่อทำทุกอย่างเสร็จแล้ว ข้อมูลก่อนสำหรับสร้างโมเดลจะได้เหมือนรูปด้านล่าง

เริ่มเขียน Code ตาม Assumption ที่ตั้งไว้ โดยกำหนดให้ p, d, q, P, D, Q มีค่าตั้งแต่ 0–2 (0,1,2) และ seasonal = 7 ใช้ itertools.product() เพื่อสร้าง pdq Components ที่เกิดขึ้นได้ทั้งหมดคือ (0,0,0), (0,0,1), …, (2,2,2)

เช่นเดียวกันกับ PDQ แต่ PDQ ยังต้องใส่ seasonal เข้าไปเป็นสมาชิกอีก 1 ตัว ดังนั้นจึงต้องเขียน List comprehension ขึ้นมาเพื่อเพิ่ม 7 เข้าไป ตัวแปร PDQS จะได้ออกมาคือ (0,0,0,7), (0,0,1,7), …, (2,2,2,7) ลองคำนวณตัวเลขง่าย ๆ สำหรับการทำ Grid-search pfq มีสมาชิก 27 ตัวซึ่งเท่ากับ PDQS ดังนั้น 27 * 27 = 729 Components จึงจำเป็นต้องใช้ Colab ในการรัน

เขียนลูปขึ้นมาเพื่อให้วนหาค่า AIC ที่ต่ำที่สุดเพื่อจะได้ค่า p,d,q และ P,D,Q ที่เหมาะที่สุดสำหรับโมเดล โดยซ้อนลูปเข้าไป 2 ชั้นเพราะต้องการทดลองทุก Combinations ที่เกิดขึ้นได้ ใช้งาน SARIMAX() โดยผ่านข้อมูลที่เราต้องการโมเดลเข้าไป ซึ่งในที่นี้คือราคาน้ำมัน กำหนด order และ seasonal_order ให้เปลี่ยนไปในแต่ละลูป กำหนด freq = ‘B’ เพราะข้อมูลราคาน้ำมันไม่ได้มีอัพเดตใหม่ทุกวัน ดังนั้นจึงต้องกำหนด Datetime ให้เป็น Business days เมื่อทุกตัวแปรพร้อมแล้ว สามารถสั่ง .fit() เพื่อเริ่มเทรนโมเดลได้เลย

อย่างที่บอกไว้ตอนต้นว่าต้องการค่า AIC ที่ต่ำที่สุด ดังนั้นก่อนเริ่มการทำงานของลูป จึงได้มีการสร้างตัวแปร aic ให้มีค่าเป็น None เพื่อรออัพเดตค่าที่เกิดจากลูปแรก เมื่อ model.fit() ลูปแรกเรียบร้อยแล้ว สามารถใช้ model.aic เพื่อดูค่า AIC ของลูปนั้น ๆ ใช้ If statement เพื่อตรวจว่าตัวแปร aic is None จากนั้นให้ Replace aic ด้วยค่า AIC จากโมเดล

Condition ที่ 2 คือถ้าตัวแปร aic ไม่เป็น None และ aic ที่เกิดขึ้นไปแล้วมีค่ามากกว่า AIC ที่ออกจากโมเดล ณ ลูปนี้ ให้อัพเดตค่าในตัวแปร aic ใหม่ด้วยค่า AIC จากโมเดลล่าสุด และเก็บค่า p,d,q และ P,D,Q,S ไว้ในตัวแปร best_pdq และ best_PDQS ตามลำดับ

สุดท้ายถ้า aic ที่มีอยู่แล้วมีค่าน้อยกว่า AIC ที่เกิดขึ้นใหม่ ให้ continue ลูปนี้ต่อไป เพราะต้องการดูทุก Combinations ที่มีโอกาสเกิดขึ้นต่อ ๆ ไป

Lowest AIC: 6520.492519031786
Best p, d, q: (1, 1, 2)
Best P, D, Q, S: (0, 1, 2, 7)

ค่า Order ทั้งสองชุดที่นำไปใช้ในโมเดลสุดท้าย

สั่งรัน SARIMAX() อีกรอบด้วยค่า best_pdq และ best_PDQS และสั่ง .fit() ตามลำดับ ถึงตอนนี้โมเดลก็พร้อมใช้งานแล้ว

Result

สามารถทดสอบกับข้อมูลอดีตได้ด้วย model.predict() เพื่อทำ Back testing เนื่องจากตัวแปร d มีค่าเท่ากับ 1 หมายถึงการใช้ First difference ทำให้ข้อมูลที่จุดแรกหายไป ดังนั้นขอให้กำหนดให้มีค่าเท่ากับข้อมูลจริงตำแหน่งแรกในอดีต กราฟด้านล่างคือการเปรียบเทียบข้อมูลจริง และผลลัพธ์ที่ได้จากโมเดล

ทดลองใช้โมเดลในการหาราคาน้ำมันในอีก 2 สัปดาห์ถัดไปด้วย Confidence interval (ความเชื่อมั่น) ที่ 90% เพื่อให้ได้ทั้งราคาที่เป็นตัวเลขเดียว และช่วงราคา Lower และ Upper โดยการ Forecast แบบ Confidence interval ให้ใช้คำสั่ง .get_forecast(steps =) จากนั้นใช้ตัวแปรเดิมด้วยคำสั่ง .conf_int(alpha =) ค่าที่ Return ออกมาคือ DataFrame ที่มีช่วงราคาอยู่ หากต้องการ Forecast เป็นเลขตัวเดียวสามารถใช้คำสั่ง .forecast(steps =)

ค่า forecast ที่ได้ออกมาจากโมเดล

ทดลองต่อข้อมูลที่เกิดขึ้นจริง กับข้อมูลที่โมเดลบอกว่ากำลังจะเกิดขึ้นในอีก 2 สัปดาห์เพื่อดูแนวโน้มที่เป็นไปได้

The oil price in the next 2 weeks is: 45.06607544102673

เพียงเท่านี้โมเดล SARIMA ถึงว่าเสร็จเป็นที่เรียบร้อย แต่ว่าโมเดลประเภทนี้สามารถเพิ่ม Assumption เข้าไปได้อีกมากมาย ซึ่งถ้าโอกาสหน้าคงมาลองเพิ่มค่าบางอย่าง แล้วดูว่าโมเดลเปลี่ยนไปอย่างไร

statisticalModel/SARIMAModel.ipynb at main · naenumtou/statisticalModel

สำหรับ Colab notebook สามารถดูได้ที่ Link ด้านบน