Category: Data analytics

readxl และ XLConnect: วิธีใช้ 2 packages สำหรับทำงานกับ Excel ในภาษา R — ตัวอย่างการทำงานกับ Daily Household Transactions

Excel เป็นเครื่องมือทำงานยอดนิยมในการทำงาน ซึ่งทำให้ในหลาย ๆ ครั้ง ข้อมูลที่เราต้องการถูกเก็บอยู่ในไฟล์ Excel (เช่น .xls และ .xlsx)

ในบทความนี้ เราจะมาทำความรู้จักกับ readxl และ XLConnect ซึ่งเป็น packages สำหรับทำงานกับ Excel ในภาษา R กัน

เราจะดูการใช้งานผ่านตัวอย่างการทำงานกับ Daily Transactions Dataset จาก Kaggle ที่ถูกเก็บในไฟล์ XLSX (”Daily Household Transactions.xlsx”):

ข้อมูลใน ”Daily Household Transactions.xlsx”

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

1️⃣ Package 1. readxl

readxl เป็น package สำหรับโหลดข้อมูลจาก Excel และมี function ที่เราจะเรียกใช้งานได้ คือ read_excel() ซึ่งต้องการ 2 arguments:

path: ชื่อไฟล์ หรือ file path
sheet: ชื่อหรือลำดับของ sheet ที่เก็บข้อมูลที่เราต้องการ

ในการเริ่มต้นใช้งาน readxl เราจะติดตั้งและโหลด package:

# Install readxl
install.packages("readxl")

# Load readxl
library(readxl)

จากนั้น เรียกใช้ read_excel() เพื่อโหลดข้อมูล:

# Import Excel data with read_excel()
all_transactions <- read_excel("Daily Household Transactions.xlsx",
                               sheet = 1)

# View the first few rows
head(all_transactions)

ผลลัพธ์:

# A tibble: 6 × 8
  Date                Mode              Category Subcategory Note  Amount `Income/Expense` Currency
  <dttm>              <chr>             <chr>    <chr>       <chr>  <dbl> <chr>            <chr>   
1 2018-09-20 12:04:08 Cash              Transpo… Train       2 Pl…     30 Expense          INR     
2 2018-09-20 12:03:15 Cash              Food     snacks      Idli…     60 Expense          INR     
3 2018-09-19 00:00:00 Saving Bank acco… subscri… Netflix     1 mo…    199 Expense          INR     
4 2018-09-17 23:41:17 Saving Bank acco… subscri… Mobile Ser… Data…     19 Expense          INR     
5 2018-09-16 17:15:08 Cash              Festiva… Ganesh Puj… Gane…    251 Expense          INR     
6 2018-09-15 06:34:17 Credit Card       subscri… Tata Sky    Perm…    200 Expense          INR

Note: read_excel() มี parametres อื่น ๆ ที่เราสามารถตั้งค่าได้ เช่น:

range: ช่วงข้อมูลที่เราต้องการโหลด
col_names: ชื่อ columns
col_types: ประเภทข้อมูลในแต่ละ column
skip: จำนวน rows ที่เราจะข้ามในการโหลดข้อมูล เช่น เรากำหนด skip = 5, read_excel() จะโหลดข้อมูลตั้งแต่ row ที่ 6 เป็นต้นไป

ดูคู่มือการใช้งาน read_excel() ทั้งหมดได้ที่ read_excel: Read xls and xlsx files.

2️⃣ Package 2. XLConnect

XLConnect เป็น package ที่ช่วยให้เราทำงานกับไฟล์ Excel จาก R ได้โดยตรง

การใช้งาน XLConnect แบ่งได้เป็น 5 กรณี ดังนี้:

โหลดและสำรวจ workbook
โหลดข้อมูลจาก sheet
จัดการ sheet
เพิ่มข้อมูลใน sheet
บันทึก workbook

เราไปดูการใช้งานในแต่ละกรณีกัน

📔 กรณีที่ 1. โหลดและสำรวจ Workbook

ในการเริ่มใช้งาน XLConnect เราจะต้องติดตั้งและโหลด package ก่อน:

# Install
install.packages("XLConnect")

# Load
library(XLConnect)

จากนั้น เราสามารถโหลด Excel เราเข้ามาใน R ได้ด้วย loadWorkbook():

# Load the workbook
workbook <- loadWorkbook("Daily Household Transactions.xlsx")

และดู sheet ทั้งหมดใน workbook ด้วย getSheets():

# List sheets
getSheets(workbook)

ผลลัพธ์:

[1] "All Transactions"

ในตัวอย่าง จะเห็นว่า Excel ของเรามี 1 sheet ได้แก่ “All Transactions”

⬇️ กรณีที่ 2. โหลดข้อมูลจาก sheet

เมื่อเห็นโครงสร้างของไฟล์ Excel แล้ว เราสามารถโหลดข้อมูลจาก sheet ที่ต้องการได้ด้วย readWorksheet():

# Get sheet data
sheet1_data <- readWorksheet(workbook,
                             sheet = "All Transactions")

# Print data
head(sheet1_data)

ผลลัพธ์:

                 Date                  Mode       Category             Subcategory                                     Note Amount Income.Expense Currency
1 2018-09-20 12:04:08                  Cash Transportation                   Train                     2 Place 5 to Place 0     30        Expense      INR
2 2018-09-20 12:03:15                  Cash           Food                  snacks              Idli medu Vada mix 2 plates     60        Expense      INR
3 2018-09-19 00:00:00 Saving Bank account 1   subscription                 Netflix                     1 month subscription    199        Expense      INR
4 2018-09-17 23:41:17 Saving Bank account 1   subscription Mobile Service Provider                        Data booster pack     19        Expense      INR
5 2018-09-16 17:15:08                  Cash      Festivals            Ganesh Pujan                              Ganesh idol    251        Expense      INR
6 2018-09-15 06:34:17           Credit Card   subscription                Tata Sky Permanent Residence - Tata Play recharge    200        Expense      INR

🖐️ กรณีที่ 3. จัดการ Sheet

ในการจัดการ sheet เราสามารถทำได้ 3 อย่าง:

สร้าง sheet ใหม่: createSheet()
เปลี่ยนชื่อ sheet: renameSheet()
ลบ sheet: removeSheet()

ยกตัวอย่างการสร้าง sheet ใหม่:

# Create new sheets
createSheet(workbook,
            name = "New")

# List sheets
getSheets(workbook)

ผลลัพธ์:

[1] "All Transactions" "New"

เปลี่ยนชื่อ sheet:

# Rename the new sheet
renameSheet(workbook,
            sheet = "New",
            newName = "Some Transactions")

# List sheets
getSheets(workbook)

ผลลัพธ์:

[1] "All Transactions"  "Some Transactions"

และลบ sheet ทิ้ง:

# Delete the new sheet
removeSheet(workbook,
            sheet = "Some Transactions")

# List sheets
getSheets(workbook)

ผลลัพธ์:

[1] "All Transactions"

➕ กรณีที่ 4. เพิ่มข้อมูลใน Sheet

เราสามารถใส่ข้อมูลลงใน sheet ได้ด้วย writeWorksheet()

ยกตัวอย่างเช่น เราต้องการเพิ่มข้อมูลสรุปค่าใช้จ่ายตามประเภทการใช้จ่าย

เราจะเริ่มจากสรุปค่าใช้จ่ายตามประเภทโดยใช้ dplyr package:

# Load dplyr
library(dplyr)

# Calculate expense by category
expense_by_cat <- sheet1_data |>
  
  # Filter for expense
  filter(Income.Expense == "Expense") |>
  
  # Group by category
  group_by(Category) |>
  
  # Calculate sum amount
  summarise(Sum = sum(Amount)) |>
    
  # Ungroup
  ungroup() |>
  
  # Sort by category
  arrange(desc(Sum))

# View the results
expense_by_cat

(Note: อ่านเกี่ยวกับการใช้ dplyr ได้ที่นี่)

ผลลัพธ์:

# A tibble: 27 × 2
   Category                  Sum
   <chr>                   <dbl>
 1 Money transfer        606529.
 2 Investment            271858 
 3 Transportation        169054.
 4 Household             161646.
 5 subscription          114588.
 6 Food                   96403.
 7 Public Provident Fund  90000 
 8 Other                  87025.
 9 Family                 78582.
10 Health                 66253.
# ℹ 17 more rows
# ℹ Use `print(n = ...)` to see more rows

จากนั้น สร้าง sheet ใหม่เพื่อเก็บข้อมูลที่เราได้มา:

# Create a new sheet
createSheet(workbook,
            name = "Expense by Category")

# Add data to "Expense by Catogory" sheet
writeWorksheet(workbook,
               data = expense_by_cat,
               sheet = "Expense by Category")

เมื่อเราเรียกดูข้อมูลจาก “Expense by Catogory” เราจะเห็นข้อมูลแบบนี้:

# Read the sheet
readWorksheet(workbook,
              sheet = "Expense by Category")

ผลลัพธ์:

                Category       Sum
1         Money transfer 606528.90
2             Investment 271858.00
3         Transportation 169053.78
4              Household 161645.58
5           subscription 114587.91
6                   Food  96403.10
7  Public Provident Fund  90000.00
8                  Other  87025.28
9                 Family  78582.20
10                Health  66252.75
11               Tourism  63608.85
12                  Gift  40168.00
13               Apparel  25373.82
14     Recurring Deposit  22000.00
15                  maid  21839.00
16                  Cook  12443.00
17                  Rent  10709.00
18             Festivals   6911.00
19               Culture   4304.36
20                Beauty   4189.00
21      Self-development   2357.00
22             Education    537.00
23              Grooming    400.00
24           Social Life    298.00
25   water (jar /tanker)    148.00
26             Documents    100.00
27      garbage disposal     67.00

💾 กรณีที่ 5. บันทึก Workbook

จนถึงจุดนี้ สิ่งที่เราแก้ไขไปจะยังไม่ถูกบันทึกลงในไฟล์ Excel ของเรา

เราต้องใช้ saveWorkbook() เพื่อบันทึกการเปลี่ยนแปลงทั้งหมดได้

อย่างในตัวอย่าง เราจะบันทึกการเปลี่ยนแปลงทั้งหมดไปที่ workbook ใหม่ชื่อ “Daily Household Transactions (Updated).xlsx”:

# Save the file
saveWorkbook(workbook,
             file = "Daily Household Transactions (Updated).xlsx")

เมื่อเราเปิด working directory ดู เราจะเห็น Excel 2 ไฟล์:

ไฟล์ต้นฉบับ
ไฟล์ที่เราสร้างใหม่

Workbook ต้นฉบับ (ซ้าย) และ workbook ที่เราสร้างใหม่ (ขวา)

💪 Summary

ในบทความนี้ เราได้เรียนรู้วิธีการทำงานกับ Excel ในภาษา R ผ่าน 2 packages ได้แก่:

Package 1. readxl:

Function	For
`read_excel()`	โหลดข้อมูลจาก Excel

Package 2. XLConnect:

Function	For
`loadWorkbook()`	โหลด workbook
`getSheets()`	ดู sheets ใน workbook
`readWorksheet()`	โหลดข้อมูลจาก sheet
`createSheet()`	สร้าง sheet
`renameSheet()`	เปลี่ยนชื่อ sheet
`removeSheet()`	ลบ sheet
`writeWorksheet()`	เพิ่มข้อมูลใน sheet
`saveWorkbook()`	บันทึก workbook

😺 GitHub

ดู code และ Excel ตัวอย่างได้ที่ GitHub

📃 References

Introduction to Importing Data in R

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

🙋 ใครที่กำลังเรียนจิตวิทยาหรือทำงานสายจิตวิทยา และเบื่อที่ต้องใช้ software ราคาแพงอย่าง SPSS และ Excel เพื่อทำข้อมูล

💪 ผมขอแนะนำ R Book for Psychologists หนังสือสอนใช้ภาษา R เพื่อการวิเคราะห์ข้อมูลทางจิตวิทยา ที่เขียนมาเพื่อนักจิตวิทยาที่ไม่เคยมีประสบการณ์เขียน code มาก่อน

ในหนังสือ เราจะปูพื้นฐานภาษา R และพาไปดูวิธีวิเคราะห์สถิติที่ใช้บ่อยกัน เช่น:

Correlation
t-tests
ANOVA
Reliability
Factor analysis

🚀 เมื่ออ่านและทำตามตัวอย่างใน R Book for Psychologists ทุกคนจะไม่ต้องพึง SPSS และ Excel ในการทำงานอีกต่อไป และสามารถวิเคราะห์ข้อมูลด้วยตัวเองได้ด้วยความมั่นใจ

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-09-18

if, for, while ใน Python: วิธีใช้ conditional statements, control flow, และ loop control ใน Python พร้อมตัวอย่าง

ในบทความนี้ เราจะมาดูวิธีใช้ code 3 ประเภท ที่จะช่วยให้ Python code สามารถตัดสินใจแทนเราได้:

Conditional statements: if, elif, else
Control flow statements: for, while
Loop control statements: continue, break, pass

ก่อนไปดูวิธีใช้ทั้ง 3 ประเภท เราจะไปทำความรู้จักกับ comparison และ logical statement ซึ่งเราจะใช้ร่วมกับ code 3 ประเภทนี้กัน

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

🧮 1. Comparison & Logical Operators

#️⃣ (1) Comparison Operators

Comparison operators เป็นเครื่องหมายใช้เปรียบเทียบข้อมูล:

Operator	Description
`==`	เท่ากับ
`!=`	ไม่เท่ากับ
`>`	มากกว่า
`>=`	มากกว่า/เท่ากับ
`<`	น้อยกว่า
`<=`	น้อยกว่า/เท่ากับ

โดย comparison operators จะให้ผลลัพธ์เป็น Boolean (True, False) กลับมา เช่น:

# True statement
10 > 5

ผลลัพธ์:

True

และ:

# False statement
10 < 5

ผลลัพธ์:

False

♟️ (2) Logical Operators

Logical operators เป็น keywords เชื่อมเงื่อนไข ช่วยให้เราประเมินหลายเงื่อนไขพร้อมกันได้:

Operator	Description
`and`	Logical AND
`or`	Logical OR
`not`	Logical NOT

โดยผลลัพธ์เป็นไปตาม truth table:

Condition 1	Operator	Condition 2	Result
`True`	`and`	`True`	`True`
`True`	`and`	`False`	`False`
`False`	`and`	`False`	`False`
`True`	`or`	`True`	`True`
`True`	`or`	`False`	`True`
`False`	`or`	`False`	`False`
–	`not`	`True`	`False`
–	`not`	`False`	`True`

ยกตัวอย่างเช่น:

# True and True
1 == 1 and 2 < 4

ผลลัพธ์:

True

🚦 2. Conditional Statements

Conditional statements ใช้กำหนดเงื่อนไขว่า code จะรันได้เมื่อไร และมีอยู่ 3 แบบ ได้แก่:

if: กำหนดเงื่อนไขแรก
elif: กำหนดเงื่อนไขเพิ่มเติม
else: รันเมื่อตรงกับเงื่อนไขอื่น ๆ ที่ไม่ใช่ if และ elif

ยกตัวอย่างเช่น เราต้องการ print ข้อความแจ้งเตือนตามสภาพอากาศ:

# Weather today
weather = "snowy"

เราสามารถทำได้แบบนี้:

# Print when sunny
if weather == "sunny":
    print("It's a sunny day. Don't forget your sunscreen!")

# Print when rainy
elif weather == "rainy":
    print("It's raining. Remember to bring an umbrella!")

# Print when other conditions
else:
    print("Likely chilly. Wear a jacket!")

ผลลัพธ์:

Likely chilly. Wear a jacket!

อธิบาย code:

if block: ประเมินว่า weather เป็น "sunny" ไหม ถ้าใช่ จะ print "It's a sunny day. Don't forget your sunscreen!"
elif block: weather เป็น "rainy" ไหม ถ้าใช่ จะ print "It's raining. Remember to bring an umbrella!"
else block: ถ้า weather เป็นค่าอื่น ๆ (เช่น "snowy") จะ print "Likely chilly. Wear a jacket!"

เพราะ weather มีค่าตรงกับ else block เราจึงได้ผลลัพธ์เป็น "Likely chilly. Wear a jacket!"

🌊 3. Control Flow Statements

Control flow statements ใช้ควบคุมลำดับการทำงานของ code และมีอยู่ 2 แบบ ได้แก่:

for: รัน code ตามจำนวนข้อมูลที่มี
while: รัน code จนกว่าเงื่อนไขจะไม่เป็นจริง

📦 (1) for Loop

ตัวอย่าง for loop เช่น print ชื่อแขกในงาน:

# Guest list
guests = ["James Bond", "John Wick", "Jack Reacher", "Jason Bourne", "Jack Ryan"]

# Print guest names
for name in guests:
    
    # Print name
    print(name)

ผลลัพธ์:

James Bond
John Wick
Jack Reacher
Jason Bourne
Jack Ryan

อธิบาย code:

guest = ...: สร้าง list เก็บรายชื่อแขก
for ...: นำชื่อแขกใน list มา print จนกว่าจะครบทุกคน

🔁 (2) while Loop

ตัวอย่าง while loop เช่น นับเลข 1 ถึง 10:

# Starting number
number = 1

# Count to 10
while number <= 10:
    
    # Print number
    print(number)
    
    # Add 1 to number
    number += 1

ผลลัพธ์:

อธิบาย code:

number = 1: กำหนดเลขเริ่มต้น
while ...: print number และปรับ number ให้สูงขึ้น 1 ค่า ทำอย่างนี้วนไปจนกว่า number จะมากกว่า 10 (number <= 10 ไม่เป็นจริง)

🚄 4. Loop Control Statements

Loop control statements ใช้ควบคุมการทำงานของ control flow และมีอยู่ 3 แบบ ได้แก่:

continue: skip ไป item ถัดไป
break: หยุดการทำงานของ control flow
pass: placeholder สำหรับใส่ code ในอนาคต

ยกตัวอย่างเช่น เรามีรายการของที่ต้องซื้อ และแต่ละอย่างมี action ไม่เหมือนกัน:

# A shopping list program
shopping_list = ["milk", "bread", "chips", "apple", "toothpaste", "chocolate"]

# Loop through the list
for item in shopping_list:

    # Skip item if chip
    if item == "chips":
        print("Chips are unhealthy. Skipping ...")
        continue
    
    # Stop the loop if toothpaste
    if item == "toothpaste":
        print("Found toothpaste, done shopping early!")
        break
    
    # Do nothing if milk
    if item == "milk" or item == "bread":
        pass
    
    # Print item
    print("Putting", item, "into the cart.")

ผลลัพธ์:

Putting milk into the cart.
Putting bread into the cart.
Chips are unhealthy. Skipping ...
Putting apple into the cart.
Found toothpaste, done shopping early!

อธิบาย code:

shopping_list = ...: สร้าง list เก็บรายการซื้อของ
if item == "chips" block: ประเมินว่า item ใช่ "chip" ไหม ถ้าใช่ ให้ print "Chips are unhealthy. Skipping ..." และข้ามไป item ถัดไป (continue)
if item == "toothpaste" block: ประเมินว่า item ใช่ "toothpaste" ไหม ถ้าใช่ ให้ print "Found toothpaste, done shopping early!" และจบ loop ทันที (break)
if item == "milk" … block: ประเมินว่า item ใช่ "milk" หรือ "bread" ไหม ถ้าใช่ ไม่ต้องทำอะไร (pass)
print ...: print ว่า กำลังเอา item ใส่ตระกร้า

สังเกตว่า:

มีแค่ milk, bread, apple ที่ “กำลังใส่ตระกร้า” เพราะ milk และ bread อยู่ใน pass block และ apple ไม่ได้อยู่ใน block ไหนเลย
Chip ถูกข้ามไป เพราะอยู่ใน continue block
Toothpaste อยู่ใน break block และมาก่อน chocolate ทำให้เราไม่เห็น chocolate เพราะ toothpaste

💪 5. Summary

ในบทความนี้ เราได้ไปดูวิธีการเขียน code ใน Python เพื่อให้ code ตัดสินใจได้:

Comparison operators:

Operator	Description
`==`	เท่ากับ
`!=`	ไม่เท่ากับ
`>`	มากกว่า
`>=`	มากกว่า/เท่ากับ
`<`	น้อยกว่า
`<=`	น้อยกว่า/เท่ากับ

Logical operators:

Operator	Description
`and`	Logical AND
`or`	Logical OR
`not`	Logical NOT

Conditional statements:

Statement	Description
`if`	กำหนดเงื่อนไขแรก
`elif`	กำหนดเงื่อนไขเพิ่มเติม
`else`	ทำ action เมื่อเข้าเงื่อนไขอื่น ๆ

Control flow statements:

Statement	Description
`for`	วนจนครบทุก item
`while`	วนจนกว่าเงื่อนไขจะเป็น `False`

Loop control statements:

Statement	Description
`continue`	Skip
`break`	Stop
`pass`	Do nothing

😺 GitHub

ดูตัวอย่าง code ทั้งหมดได้ที่ GitHub

📃 References

2025-09-11

4 วิธีการใช้ getSymbols() เพื่อโหลดข้อมูลการเงินในภาษา R — ตัวอย่างการโหลดข้อมูลหุ้น Apple

ในบทความนี้ เราจะมาทำความรู้จักกับ getSymbols() จาก quantmod package ซึ่งใช้โหลดข้อมูลทางการเงิน เช่น ข้อมูลหุ้น และข้อมูลทางเศรษฐกิจ กัน

โดยเราจะไปดู 4 วิธีการใช้งาน ได้แก่:

Basics: การใช้งาน getSymbols() เบื้องต้น
Advanced: การใช้งานขั้นสูง
View columns: การดู column ข้อมูล
Visualise data: การสร้างกราฟข้อมูล

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

💻 Install & Import

ในการใช้งาน getSymbols() เราจะต้องเริ่มจากการติดตั้งและโหลด quantmod package ก่อน แบบนี้:

# Install the packages
install.packages("quantmod")

# Load the packages
library(quantmod)

1️⃣ Basics: 5 Parametres to Know

ในการใช้งานพื้นฐาน getSymbols() มี 5 parametres หลักที่เราควรรู้ ได้แก่:

Symbols: ตัวย่อชื่อข้อมูล
src: แหล่งข้อมูล
auto.assign: ให้โหลด (TRUE) หรือแสดงข้อมูล (FALSE)
env: environment สำหรับโหลดข้อมูล
from และ to: ช่วงเวลาของข้อมูล

เราไปดูวิธีใช้งานแต่ละ parametre กัน

💲 Symbols

ในขั้นแรกของการโหลดข้อมูลการเงินด้วย getSymbols() ให้เราระบุชื่อข้อมูลที่เราต้องการ เช่น:

ข้อมูล	ชื่อข้อมูล
Apple	AAPL
Google	GOOG
Microsoft	MSFT
NVIDIA	NVDA

ในบทความนี้ เราจะลองโหลดข้อมูลหุ้น Apple กัน ซึ่งเราสามารถกำหนด Symbols ได้แบบนี้:

# Import Apple data
getSymbols("AAPL")

เมื่อรันแล้ว เราจะได้ตัวแปรประเภท xts (time series object) ชื่อ AAPL มา ซึ่งเราดูข้อมูลในตัวแปรนี้ได้ด้วย head():

# Print result
head(AAPL)

ผลลัพธ์:

           AAPL.Open AAPL.High AAPL.Low AAPL.Close AAPL.Volume AAPL.Adjusted
2007-01-03  3.081786  3.092143 2.925000   2.992857  1238319600      2.518541
2007-01-04  3.001786  3.069643 2.993571   3.059286   847260400      2.574443
2007-01-05  3.063214  3.078571 3.014286   3.037500   834741600      2.556108
2007-01-08  3.070000  3.090357 3.045714   3.052500   797106800      2.568732
2007-01-09  3.087500  3.320714 3.041071   3.306071  3349298400      2.782115
2007-01-10  3.383929  3.492857 3.337500   3.464286  2952880000      2.915257

จะเห็นได้ว่า นอกจากวันที่แล้ว ข้อมูลเรายังประกอบไปด้วย 6 columns ดังนี้:

Open: ราคาเปิด
High: ราคาสูงสุด
Low: ราคาต่ำสุด
Close: ราคาปิด
Volume: ปริมาณการซื้อขาย
Adjusted: ราคาปรับปรุง

ในกรณีที่เราต้องการโหลดข้อมูลหุ้นหลายตัวพร้อมกัน เราสามารถใช้ character vector ช่วยได้แบบนี้:

# Load multiple instruments
getSymbols(c("AAPL", "GOOGL", "MSFT", "NVDA"))

🏦 src

src ใช้กำหนดแหล่งข้อมูลที่ getSymbols() จะไปดึงข้อมูลมา

โดย default, src ถูกตั้งไว้ให้ดึงข้อมูลจาก Yahoo! Finance ("yahoo")

แสดงว่า เราสามารถเขียน code เพื่อดึงข้อมูล Apple จาก Yahoo! Finance ได้ทั้งแบบนี้:

# Import Apple data
getSymbols("AAPL")

และแบบนี้:

# Import Apple data
getSymbols("AAPL", src = "yahoo")

ทั้งนี้ getSymbols() มีแหล่งข้อมูลอื่น ๆ ให้เราเลือกใช้งานได้ เช่น:

Federal Reserve Economic Data ("FRED")
Oanda ("oanda")

นอกจากแหล่งข้อมูลออนไลน์แล้ว เรายังสามารถโหลดข้อมูลแบบออฟไลน์ได้ เช่น:

CSV ("csv")
RData ("RData")

ยกตัวอย่างเช่น โหลดข้อมูล Apple จาก CSV:

# Load csv data
getSymbols("AAPL", src = "csv")

Note: การโหลดข้อมูล CSV ด้วย getSymbols() มีเงื่อนไข 4 อย่าง ได้แก่:

มีชื่อไฟล์เป็นชื่อหุ้น (เช่น AAPL.csv)
มี header
Column แรกจะต้องเป็น datetime (เช่น 2025-05-25)
Column ที่เหลือควรมีชื่อและข้อมูลดังนี้:
1. Open: ข้อมูลราคาเปิด
2. High: ข้อมูลราคาสูงสุด
3. Low: ข้อมูลราคาต่ำสุด
4. Close: ข้อมูลราคาปิด
5. Volume: ข้อมูลปริมาณการซื้อขาย
6. Adjusted: ข้อมูลราคาปรับปรุง

🏧 auto.assign

auto.assign ใช้กำหนดว่า getSymbols() จะโหลดข้อมูลมาไว้ใน global environment หรือแสดงข้อมูลที่โหลดมาได้

โดย default, auto.assign เป็น TRUE ซึ่งทำให้เราได้ตัวแปรที่เก็บข้อมูลการเงินมาไว้ใน global environment ของเรา โดยไม่ต้องกำหนดตัวแปรเอง

ทั้งนี้ ในกรณีที่เราต้องการกำหนดตัวแปรเอง ให้เราเปลี่ยน auto.assign เป็น FALSE แบบนี้:

# Set auto.assign to FALSE to assign to custom variable
apple_data <- getSymbols("AAPL", auto.assign = FALSE)

ถ้าเรากำหนดให้ auto.assign = FALSE โดยไม่จัดเก็บไว้ในตัวแปร getSymbols() จะแสดงข้อมูลใน console ของเรา:

# Set auto.assign to FALSE without variable assignment
getSymbols("AAPL", auto.assign = FALSE)

ผลลัพธ์:

            AAPL.Open  AAPL.High   AAPL.Low AAPL.Close AAPL.Volume AAPL.Adjusted
2007-01-03   3.081786   3.092143   2.925000   2.992857  1238319600      2.518541
2007-01-04   3.001786   3.069643   2.993571   3.059286   847260400      2.574442
2007-01-05   3.063214   3.078571   3.014286   3.037500   834741600      2.556108
2007-01-08   3.070000   3.090357   3.045714   3.052500   797106800      2.568732
2007-01-09   3.087500   3.320714   3.041071   3.306071  3349298400      2.782115
2007-01-10   3.383929   3.492857   3.337500   3.464286  2952880000      2.915257
2007-01-11   3.426429   3.456429   3.396429   3.421429  1440252800      2.879192
2007-01-12   3.378214   3.395000   3.329643   3.379286  1312690400      2.843727
2007-01-16   3.417143   3.473214   3.408929   3.467857  1244076400      2.918262
2007-01-17   3.484286   3.485714   3.386429   3.391071  1646260000      2.853645
       ...                                                                      
2025-06-05 203.500000 204.750000 200.149994 200.630005    55126100    200.630005
2025-06-06 203.000000 205.699997 202.050003 203.919998    46607700    203.919998
2025-06-09 204.389999 206.000000 200.020004 201.449997    72862600    201.449997
2025-06-10 200.600006 204.350006 200.570007 202.669998    54672600    202.669998
2025-06-11 203.500000 204.500000 198.410004 198.779999    60989900    198.779999
2025-06-12 199.080002 199.679993 197.360001 199.199997    43904600    199.199997
2025-06-13 199.729996 200.369995 195.699997 196.449997    51447300    196.449997
2025-06-16 197.300003 198.690002 196.559998 198.419998    43020700    198.419998
2025-06-17 197.199997 198.389999 195.210007 195.639999    38856200    195.639999
2025-06-18 195.940002 197.570007 195.070007 196.580002    45350400    196.580002

🌲 env

env ใช้กำหนด environment ที่ใช้เก็บข้อมูล ซึ่งโดย default, getSymbols() จะโหลดข้อมูลไว้ใน global environment

ในการใช้ env กำหนด environment ที่ต้องการ เราจะต้องเริ่มจากสร้าง environment ขึ้นมาก่อนด้วย new.env():

# Create a new environment
my_env <- new.env()

จากนั้น โหลดข้อมูลเข้าไปใน environment ใหม่:

# Load data into the environment
getSymbols("AAPL", env = my_env)

เราสามารถดูตัวแปรที่เก็บไว้ใน environment ได้ด้วย ls():

# List all variables in environment
ls(envir = my_env)

และดูข้อมูลได้ด้วย $ เช่น:

# Show Apple data
head(my_env$AAPL)

ผลลัพธ์:

           AAPL.Open AAPL.High AAPL.Low AAPL.Close AAPL.Volume AAPL.Adjusted
2007-01-03  3.081786  3.092143 2.925000   2.992857  1238319600      2.518541
2007-01-04  3.001786  3.069643 2.993571   3.059286   847260400      2.574442
2007-01-05  3.063214  3.078571 3.014286   3.037500   834741600      2.556109
2007-01-08  3.070000  3.090357 3.045714   3.052500   797106800      2.568732
2007-01-09  3.087500  3.320714 3.041071   3.306071  3349298400      2.782115
2007-01-10  3.383929  3.492857 3.337500   3.464286  2952880000      2.915256

📆 from, to

from ใช้กำหนดวันแรกของข้อมูล และ to กำหนดวันสุดท้ายของข้อมูล

ยกตัวอย่างเช่น เราต้องการโหลดข้อมูล Apple ในเดือนพฤษภาคม 2025:

# Load data for May 2025
apple_data_2025_05 = getSymbols("AAPL",
                                auto.assign = FALSE,
                                from = "2025-05-01",
                                to = "2025-05-31")

# Print results
print("First three records:")
head(apple_data_2025_05, n = 3)
print("------------------------------------------------------------------------------")
print("Last three records:")
tail(apple_data_2025_05, n = 3)

ผลลัพธ์:

First three records:
           AAPL.Open AAPL.High AAPL.Low AAPL.Close AAPL.Volume AAPL.Adjusted
2025-05-01    209.08    214.56   208.90     213.32    57365700      213.0406
2025-05-02    206.09    206.99   202.16     205.35   101010600      205.0811
2025-05-05    203.10    204.10   198.21     198.89    69018500      198.6295
------------------------------------------------------------------------------
Last three records:
           AAPL.Open AAPL.High AAPL.Low AAPL.Close AAPL.Volume AAPL.Adjusted
2025-05-28    200.59    202.73   199.90     200.42    45339700        200.42
2025-05-29    203.58    203.81   198.51     199.95    51396800        199.95
2025-05-30    199.37    201.96   196.78     200.85    70819900        200.85

Note: ตลาดหุ้นวันสุดท้ายของเดือนพฤษภาคม 2025 คือ 30 พฤษภาคม ทำให้ข้อมูลสิ้นสุด ณ วันที่ 30

2️⃣ Advanced: Set Defaults

จะเห็นว่า getSymbols() การตั้งค่าที่หลากหลาย

ทั้งนี้ quantmod มี 6 functions ที่ใช้ร่วมกับ getSymbols() เพื่อช่วยลดขั้นตอนในการตั้งค่าต่าง ๆ ได้แก่:

getDefaults() และ setDefaults()
getSymbolLookup() และ setSymbolLookup()
saveSymbolLookup() และ loadSymbolLookup()

เราไปดูวิธีการใช้งานทั้ง 6 functions กัน

🌍 getDefaults(), setDefaults()

getDefaults() ใช้สำหรับดูค่า default ของ getSymbols()

ส่วน setDefaults() ใช้สำหรับกำหนดค่า default

ยกตัวอย่างเช่น เราต้องการกำหนด:

src จาก "yahoo" เป็น "FRED"
auto.assign จาก TRUE เป็น FALSE

เราสามารถทำได้แบบนี้:

# Get defaults before changing
print("Defaults (before):")
getDefaults(getSymbols)

# Set defaults
setDefaults(getSymbols,
            src = "FRED",
            auto.assign = FALSE)
 
# Get defaults after changing
print("Defaults (after):")
getDefaults(getSymbols)

ผลลัพธ์:

Defaults (before):
NULL

Defaults (after):
$src
[1] "'FRED'"

$auto.assign
[1] FALSE

ตอนนี้ ถ้าเราเรียกใช้ getSymbols() โดยไม่กำหนด src และ auto.assign ทั้งสอง arguments นี้จะเป็น "FRED" และ FALSE ตามลำดับ

ถ้าเราต้องการ reset ค่า default ให้เราใส่ NULL ใน setDefaults():

# Reset defaults
setDefaults(getSymbols,
            src = NULL,
            auto.assign = NULL)

# Check defaults after resetting
print("Defaults (reset):")
getDefaults(getSymbols)

ผลลัพธ์:

Defaults (reset):
NULL

💹 setSymbolLookup(), getSymbolLookup()

getSymbolLookup() และ setSymbolLookup() ทำงานเหมือนกับ getDefaults() และ setDefaults() แต่เป็นการตั้งค่า default สำหรับข้อมูลแต่ละตัว (แทนระดับ global) ภายใน session เท่านั้น

ยกตัวอย่างเช่น เราต้องการกำหนดให้ดึงข้อมูลหุ้น Google มาจาก Google Finance:

# Set default for Google
setSymbolLookup(GOOG = list(src = "google"))

# Get new defaults
getSymbolLookup()

ผลลัพธ์:

$GOOG
$GOOG$src
[1] "google"

Note: Google Finance หยุดให้ข้อมูลกับ quantmod เมื่อปี 2018 ทำให้ตอนนี้ เราไม่สามารถใช้ src = "google" ได้อีก

เราสามารถ reset ค่า default ได้โดยใส่ NULL เหมือนเดิม:

# Reset defaults
setSymbolLookup(NULL)

# Get defaults after resetting
getSymbolLookup()

ผลลัพธ์:

named list()

💾 saveSymbolLookup(), loadSymbolLookup()

ค่าที่เราใช้ setSymbolLookup() กำหนดจะถูก reset ทุกครั้งที่เราปิด session ไป

ถ้าเราต้องการเก็บค่าเพื่อไปใช้ใน session อื่น เราสามารถใช้ saveSymbolLookup() และ loadSymbolLookup() ช่วยได้:

saveSymbolLookup() บันทึกค่าเก็บไว้ในไฟล์ RDS
loadSymbolLookup() โหลดค่าที่เก็บไว้ใน RDS

ตัวอย่าง code:

# Save defaults
saveSymbolLookup(file = "symbols.rds")

# Load defaults
loadSymbolLookup(file = "symbols.rds")

3️⃣ View Columns

หลังจากเราโหลดข้อมูลมาแล้ว เราสามารถดูข้อมูลแต่ละ column ได้ด้วย 7 functions หลัก ดังนี้:

No.	Function	For
1	`Op()`	ราคาเปิด
2	`Hi()`	ราคาสูงสุด
3	`Lo()`	ราคาต่ำสุด
4	`Cl()`	ราคาปิด
5	`Ad()`	ราคาปรับปรุง
6	`Vo()`	ปริมาณการซื้อขาย
7	`OHLC()`	ราคาเปิด, สูงสุด, ต่ำสุด, และราคาปิด

ตัวอย่างเช่น ดูราคาเปิด:

# Get opening price
Op(AAPL)

ผลลัพธ์:

            AAPL.Open
2007-01-03   3.081786
2007-01-04   3.001786
2007-01-05   3.063214
2007-01-08   3.070000
2007-01-09   3.087500
2007-01-10   3.383929
2007-01-11   3.426429
2007-01-12   3.378214
2007-01-16   3.417143
2007-01-17   3.484286
       ...           
2025-06-05 203.500000
2025-06-06 203.000000
2025-06-09 204.389999
2025-06-10 200.600006
2025-06-11 203.500000
2025-06-12 199.080002
2025-06-13 199.729996
2025-06-16 197.300003
2025-06-17 197.199997
2025-06-18 195.940002

ดูราคาสูงสุด:

# Get highest price
Hi(AAPL)

ผลลัพธ์:

            AAPL.High
2007-01-03   3.092143
2007-01-04   3.069643
2007-01-05   3.078571
2007-01-08   3.090357
2007-01-09   3.320714
2007-01-10   3.492857
2007-01-11   3.456429
2007-01-12   3.395000
2007-01-16   3.473214
2007-01-17   3.485714
       ...           
2025-06-05 204.750000
2025-06-06 205.699997
2025-06-09 206.000000
2025-06-10 204.350006
2025-06-11 204.500000
2025-06-12 199.679993
2025-06-13 200.369995
2025-06-16 198.690002
2025-06-17 198.389999
2025-06-18 197.570007

ดูราคาต่ำสุด:

# Get lowest price
Lo(AAPL)

ผลลัพธ์:

             AAPL.Low
2007-01-03   2.925000
2007-01-04   2.993571
2007-01-05   3.014286
2007-01-08   3.045714
2007-01-09   3.041071
2007-01-10   3.337500
2007-01-11   3.396429
2007-01-12   3.329643
2007-01-16   3.408929
2007-01-17   3.386429
       ...           
2025-06-05 200.149994
2025-06-06 202.050003
2025-06-09 200.020004
2025-06-10 200.570007
2025-06-11 198.410004
2025-06-12 197.360001
2025-06-13 195.699997
2025-06-16 196.559998
2025-06-17 195.210007
2025-06-18 195.070007

ดูราคาปิด:

# Get closing price
Cl(AAPL)

ผลลัพธ์:

           AAPL.Close
2007-01-03   2.992857
2007-01-04   3.059286
2007-01-05   3.037500
2007-01-08   3.052500
2007-01-09   3.306071
2007-01-10   3.464286
2007-01-11   3.421429
2007-01-12   3.379286
2007-01-16   3.467857
2007-01-17   3.391071
       ...           
2025-06-05 200.630005
2025-06-06 203.919998
2025-06-09 201.449997
2025-06-10 202.669998
2025-06-11 198.779999
2025-06-12 199.199997
2025-06-13 196.449997
2025-06-16 198.419998
2025-06-17 195.639999
2025-06-18 196.580002

ดูราคาปรับปรุง:

# Get adjusted price
Ad(AAPL)

ผลลัพธ์:

           AAPL.Adjusted
2007-01-03      2.518541
2007-01-04      2.574442
2007-01-05      2.556109
2007-01-08      2.568732
2007-01-09      2.782116
2007-01-10      2.915256
2007-01-11      2.879192
2007-01-12      2.843727
2007-01-16      2.918261
2007-01-17      2.853645
       ...              
2025-06-05    200.630005
2025-06-06    203.919998
2025-06-09    201.449997
2025-06-10    202.669998
2025-06-11    198.779999
2025-06-12    199.199997
2025-06-13    196.449997
2025-06-16    198.419998
2025-06-17    195.639999
2025-06-18    196.580002

ดูปริมาณการซื้อขาย:

# Get volume
Vo(AAPL)

ผลลัพธ์:

           AAPL.Volume
2007-01-03  1238319600
2007-01-04   847260400
2007-01-05   834741600
2007-01-08   797106800
2007-01-09  3349298400
2007-01-10  2952880000
2007-01-11  1440252800
2007-01-12  1312690400
2007-01-16  1244076400
2007-01-17  1646260000
       ...            
2025-06-05    55126100
2025-06-06    46607700
2025-06-09    72862600
2025-06-10    54672600
2025-06-11    60989900
2025-06-12    43904600
2025-06-13    51447300
2025-06-16    43020700
2025-06-17    38856200
2025-06-18    45350400

ดูราคาเปิด, สูงสุด, ต่ำสุด, และราคาปิด:

# Get all price
OHLC(AAPL)

ผลลัพธ์:

            AAPL.Open  AAPL.High   AAPL.Low AAPL.Close
2007-01-03   3.081786   3.092143   2.925000   2.992857
2007-01-04   3.001786   3.069643   2.993571   3.059286
2007-01-05   3.063214   3.078571   3.014286   3.037500
2007-01-08   3.070000   3.090357   3.045714   3.052500
2007-01-09   3.087500   3.320714   3.041071   3.306071
2007-01-10   3.383929   3.492857   3.337500   3.464286
2007-01-11   3.426429   3.456429   3.396429   3.421429
2007-01-12   3.378214   3.395000   3.329643   3.379286
2007-01-16   3.417143   3.473214   3.408929   3.467857
2007-01-17   3.484286   3.485714   3.386429   3.391071
       ...                                            
2025-06-05 203.500000 204.750000 200.149994 200.630005
2025-06-06 203.000000 205.699997 202.050003 203.919998
2025-06-09 204.389999 206.000000 200.020004 201.449997
2025-06-10 200.600006 204.350006 200.570007 202.669998
2025-06-11 203.500000 204.500000 198.410004 198.779999
2025-06-12 199.080002 199.679993 197.360001 199.199997
2025-06-13 199.729996 200.369995 195.699997 196.449997
2025-06-16 197.300003 198.690002 196.559998 198.419998
2025-06-17 197.199997 198.389999 195.210007 195.639999
2025-06-18 195.940002 197.570007 195.070007 196.580002

4️⃣ Visualise Data

สุดท้าย เราสามารถสร้างกราฟเพื่อสำรวจข้อมูล โดยใช้ 2 functions ได้แก่:

autoplot() จาก ggplot2 package
chartSeries() จาก quantmod package

ยกตัวอย่างเช่น สำรวจราคาปิดของ Apple:

# Import ggplots
library(ggplot2)

# Visualise with autoplot()
autoplot(Cl(AAPL),
         ts.colour = "darkgreen") +

  # Add text
  labs(title = "AAPL Closing Price (Jan 2007 – Jun 2025)",
       x = "Time",
       y = "Price (USD)")

ผลลัพธ์:

# Visualise with chartSeries()
chartSeries(Cl(AAPL))

ผลลัพธ์:

💪 Summary

ในบทความนี้ เราได้ดูวิธีการทำงานกับ getSymbols() เพื่อโหลดข้อมูลจากการเงินจากแหล่งต่าง ๆ กัน

ตอนนี้ เรารู้จักกับ 5 parametres หลักของ getSymbols():

Symbols
src
auto.assign
env
from และ to

วิธีตั้งค่า default ด้วย 3 คู่ functions:

getDefaults() และ setDefaults()
getSymbolLookup() และ setSymbolLookup()
saveSymbolLookup() และ loadSymbolLookup()

วิธีดูข้อมูลด้วย 7 functions:

Op()
Hi()
Lo()
Cl()
Ad()
Vo()
OHLC()

และสุดท้าย วิธีสร้างกราฟด้วย 2 functions:

autoplot()
chartSeries()

😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-08-28

Regular Expressions: เครื่องมือถอดรหัส เพื่อทำงานกับ Text ใน Google Sheets — แนะนำวิธีใช้สูตร REGEXMATCH, REGEXEXTRACT, REGEXREPLACE

Regular expression หรือ การเขียน pattern เพื่อทำงานกับ text เป็นเหมือนกับเครื่องถอดรหัสที่ช่วยถอดรูปแบบของข้อความ และช่วยให้เราทำงานกับข้อมูลได้ง่ายขึ้น

ยกตัวอย่างเช่น เรามีข้อมูลลูกค้าที่เก็บค่าไว้ใน cell เดียวกัน:

แต่เราต้องการแยกข้อมูลเป็น columns เพื่อความสะดวกในการใช้งานต่อ เช่น:

ID
Name
Address

แทนที่เรา copy-paste ข้อมูลลงทีละ column ด้วยตัวเอง เราสามารถใช้สูตร regular expression ของ Google Sheets ช่วยได้:

ในบทความนี้ เราจะทำความรู้จักกับ:

3 สูตร Google Sheets ที่ใช้ทำงานกับ regular expression
การเขียน regular expression เพื่อใช้งาน 3 สูตรนี้

ถ้าพร้อมแล้วมาเริ่มกันเลย

👨‍🔬 Intro to REGEX Formulas

Google Sheets มี 3 สูตรสำหรับใช้ regular expression ซึ่งได้แก่:

No.	Formula	For
1	`REGEXMATCH()`	เช็กหาคำที่ตรงกับ regular expression
2	`REGEXEXTRACT()`	ดึงคำที่ตรงกับ regular expression ออกมา
3	`REGEXREPLACE()`	แทนที่คำที่ตรงกับ regular expression

โดยมีการเขียนดังนี้:

No.	Formula	Syntax
1	`REGEXMATCH()`	`REGEXMATCH(text, regular_expression)`
2	`REGEXEXTRACT()`	`REGEXEXTRACT(text, regular_expression)`
3	`REGEXREPLACE()`	`REGEXREPLACE(text, regular_expression, replacement)`

text คือ ข้อมูล text ที่เราจะใช้ทำงาน
regular_expression คือ regular expression ที่เราใช้กำหนดเงื่อนไข
replacement คือ คำที่จะแทนที่คำที่ตรงกับ regular expression

ก่อนไปดูการใช้งานจริง เรามาดูหลักการเขียน regular expression เพื่อให้สูตรเหล่านี้ทำงานได้กัน

⌨️ Regular Expression Basics

การเขียน regular expression แบ่งออกได้เป็น 4 กลุ่มหลัก ได้แก่:

Special characters
Character classes
Quantifiers
Others

🔣 (1) Special Characters

Regular expression ในกลุ่มนี้ เป็นอักขระพิเศษ โดยมี usage ต่างกันดังนี้:

No.	Character	For	Example
1	`^…`	เริ่มต้นด้วย …	`^J` หมายถึง คำที่ขึ้นต้นด้วย J
2	`…$`	ลงท้ายด้วย …	`s$` หมายถึง คำที่ลงท้ายด้วย s
3	`.`	แทนอักขระใด ๆ 1 ตัว (ยกเว้นการขึ้นบรรทัดใหม่)	`c.t` ใช้หาคำเช่น cat, cmt, cPt, c7t
5	`\\…`	ใช้บอกว่า … ไม่ใช่ regular expression	`\\.` หมายถึง ให้หาคำที่มี . แทนคำที่มีอักขระใด ๆ

🔣 (2) Character Classes

Regular expression ในกลุ่มนี้ ใช้แทนตัวอักษรและตัวเลข:

No.	Character	For	Example
1	`[…]`	จับคู่ตัวอักษร ใน `[]`	`[ABC]` หมายถึง จับคู่คำที่มี A, B, หรือ C
2	`[A-Z]`	จับคู่ตัวอักษร A-Z พิมพ์ใหญ่	จะได้คำ เช่น ANT, FGH, QPD, …
3	`[a-z]`	จับคู่ตัวอักษร a-z พิมพ์เล็ก	จะได้คำ เช่น ant, fgh, qpd, …
4	`[A-Za-z]`	จับคู่ตัวอักษร A-Z พิมพ์เล็กหรือใหญ่ก็ได้	จะได้คำ เช่น aNt, Fgh, qpD, …
5	`[A-z]`	จับคู่ตัวอักษร a-z พิมพ์เล็กหรือใหญ่ก็ได้	จะได้คำ เช่น aNt, Fgh, qpD, …
6	`[0-9]`	จับคู่ตัวตัวเลข 0-9	จะได้คำ เช่น 012, 438, 507, …

🔣 (3) Quantifiers

Regular expression ในกลุ่มนี้ ใช้แทนจำนวนตัวอักขระที่ต้องการค้นหา:

No.	Character	For	Example
1	`…*`	จับคู่อักขระ … จำนวนตั้งแต่ 0 ขึ้นไป	`A*` จับคู่คำ เช่น A, An, Ant, ANts, …
2	`…+`	จับคู่อักขระ … จำนวนตั้งแต่ 1 ขึ้นไป	`A+` จับคู่คำ เช่น An, Ant, ANts, …
3	`…?`	จับคู่อักขระ … จำนวน 0 หรือ 1	`A?` จับคู่คำ เช่น A, An, AT, …
5	`…{…}`	จับคู่อักขระ … ตามจำนวนใน `{}`	`.{5}` หมายถึง จับคู่อักขระใด ๆ ที่มีความยาว 5 ตัว
4	`…{…, …}`	จับคู่อักขระ … ตามช่วงใน `{}`	`.{3, 5}` หมายถึง จับคู่อักขระใด ๆ ที่มีความยาว 3 ถึง 5 ตัว

🔣 (4) Others

Regular expression อื่น ๆ นอกเหนือ 3 กลุ่มแรก:

No.	Character	For	Example
1	`\\w`	จับคู่ A-Z ทั้งพิมพ์เล็กและใหญ่ หรือตัวเลขก็ได้	–
2	`\\d`	จับคู่ตัวตัวเลข 0-9	–
3	`\\b`	จับคู่ word boundary	–
4	`\\n`	จับคู่ new line	–
5	`\\s`	จับคู่ blank space	–
6	`(…)`	จับกลุ่ม regular expression	–

📃 Learn More

เราสามารถ regular expression อื่น ๆ ได้ที่ Syntax for Regular Expressions จาก Google

💪 Put It All Together

ตอนนี้ เรารู้จัก 3 สูตร REGEX และการเขียน regular expression แล้ว

เรามาดูตัวอย่างการใช้งานสูตร REGEX ในการทำงานกัน

ในตัวอย่าง เรามีข้อมูลพนักงาน 3 columns ได้แก่:

Name
Email
Comment

Note: สามารถดูต้นฉบับได้ที่ Google Sheets

🗝️ (1) REGEXMATCH

สมมุติว่า เราต้องการเช็กว่า ใน comment มีคำว่า “comment” ไหม

เราสามารถใช้ REGEXMATCH ได้ดังนี้:

=REGEXMATCH(C3, "[Cc]omment")

ผลลัพธ์:

จะเห็นว่า เราได้ค่า TRUE และ FALSE กลับมา:

TRUE แสดงว่า text มีคำว่า “comment”
FALSE แสดงว่า text ไม่มีคำว่า “comment”

🗝️ (2) REGEXEXTRACT

เราต้องการดึงนามสกุลออกจากชื่อ

เราสามารถใช้ REGEXEXTRACT ช่วยได้:

=REGEXEXTRACT(A3, "\s([A-z]*)")

ผลลัพธ์:

🗝️ (3) REGEXREPLACE

เราต้องการเปลี่ยนอีเมลของทุกคนให้มี domain (เช่น @example.com และ @adams-mail.net) เป็น @email.com

เราสามารถใช้ REGEXREPLACE ช่วยได้:

=REGEXREPLACE(B3, "@[A-z.-]+\.[A-z]{2,}", "@email.com")

ผลลัพธ์:

📚 References

2025-08-14

7 วิธีทำงานกับ time series data ในภาษา R: Base R, lubridate, และ zoo packages — ตัวอย่างการวิเคราะห์ราคา Bitcoin

ในบทความนี้ เราจะมาเรียนรู้ 7 วิธีการทำงานกับ time series data หรือข้อมูลที่จัดเรียงด้วยเวลา ในภาษา R ผ่านการทำงานกับ Bitcoin Historical Data ซึ่งมีข้อมูล Bitcoin ในช่วงปี ค.ศ. 2010–2024 กัน:

Converting: การแปลงข้อมูลเป็น datetime และ time series
Missing values: การจัดการ missing values ใน time series data
Plotting: การสร้างกราฟ time series data
Subsetting: การเลือกข้อมูลจาก time series data
Aggregating: การสรุป time series data
Rolling window: การทำ rolling window
Expanding window: การทำ expanding window

ถ้าพร้อมแล้วไปเริ่มกันเลย

🏁 Getting Started

ในการเริ่มทำงานกับ time series ในภาษา R เราจะเริ่มต้นจาก 2 อย่าง ได้แก่:

Install and load packages
Load dataset

1️⃣ Install & Load Packages

ในภาษา R เรามี 3 packages ที่ใช้ทำงานกับ time series data บ่อย ๆ ได้แก่:

Base R packages ที่มาพร้อมกับภาษา R
lubridate: ใช้แปลงข้อมูลและจัด format ของ date-time data
zoo: ใช้ทำงานกับ time series data

ในการเริ่มต้นใช้งาน เราจะติดตั้ง packages เหล่านี้กัน:

# Install packages
install.packages("lubridate"
install.packages("zoo")
install.packages("ggplot2") # for data visualisation

จากนั้น เรียกใช้งาน packages:

# Load packages
library(lubridate)
library(zoo)
library(ggplot2) # for data visualisation

Note: เราเพิ่ม ggplot2 เข้ามาด้วย เพื่อช่วยในการสร้างกราฟในภายหลัง

2️⃣ Load Dataset

หลังติดตั้งและเรียกใช้งาน packages เราจะโหลด Bitcoin Historical Data (ดาวน์โหลดไฟล์ได้ที่ Kaggle) เข้ามาใน environment:

# Load the dataset
btc <- read.csv("btc_hist_2010-2024.csv")

จากนั้น เราสามารถสำรวจข้อมูลได้ด้วย head() และ str():

## View the first 6 rows
head(btc)

## View the structure
str(btc)

ผลลัพธ์ของ head():

          Date    Price     Open     High      Low   Vol. Change..
1 Feb 09, 2024 47,545.4 45,293.3 47,710.2 45,254.2 86.85K    4.97%
2 Feb 08, 2024 45,293.3 44,346.2 45,579.2 44,336.4 66.38K    2.15%
3 Feb 07, 2024 44,339.8 43,088.4 44,367.9 42,783.5 48.57K    2.91%
4 Feb 06, 2024 43,087.7 42,697.6 43,375.5 42,566.8 33.32K    0.91%
5 Feb 05, 2024 42,697.2 42,581.4 43,532.2 42,272.5 39.26K    0.27%
6 Feb 04, 2024 42,581.4 43,006.2 43,113.2 42,379.4 20.33K   -0.99%

ผลลัพธ์ของ str():

'data.frame':	4955 obs. of  7 variables:
 $ Date    : chr  "Feb 09, 2024" "Feb 08, 2024" "Feb 07, 2024" "Feb 06, 2024" ...
 $ Price   : chr  "47,545.4" "45,293.3" "44,339.8" "43,087.7" ...
 $ Open    : chr  "45,293.3" "44,346.2" "43,088.4" "42,697.6" ...
 $ High    : chr  "47,710.2" "45,579.2" "44,367.9" "43,375.5" ...
 $ Low     : chr  "45,254.2" "44,336.4" "42,783.5" "42,566.8" ...
 $ Vol.    : chr  "86.85K" "66.38K" "48.57K" "33.32K" ...
 $ Change..: chr  "4.97%" "2.15%" "2.91%" "0.91%" ...

จากผลลัพธ์ เราจะเห็นว่า Bitcoin Historical Data มี 7 columns ด้วยกัน ได้แก่:

No.	Column	Description
1	`Date`	วันที่
2	`Close`	ราคาปิด (USD)
3	`Open`	ราคาเปิด (USD)
4	`High`	ราคาสูงสุด (USD)
5	`Low`	ราคาต่ำสุด (USD)
6	`Vol.`	ปริมาณการซื้อขาย
7	`Change..`	การเปลี่ยนแปลงราคาจากวันก่อน (คิดเป็น %)

นอกจากนี้ เราจะเห็นได้ว่า ทุก column เป็น character ซึ่งเราจะต้องแปลงให้เป็นประเภทข้อมูลที่ถูกต้องก่อน ซึ่งเราสามารถทำได้ดังนี้:

สร้าง copy ของ btc
เปลี่ยนชื่อ column (เพื่อให้ง่ายต่อการทำงาน)
เปลี่ยน format ข้อมูลของ price และ volume
แปลงประเภทข้อมูลจาก character เป็น numeric

# Clean the dataset

## Make a copy of the dataset
btc_cleaned <- btc

## Rename columns

### Create new column names
new_col_names <- c("date",
                   "close", "open",
                   "high", "low",
                   "volume", "change")

### Rename
colnames(btc_cleaned) <- new_col_names

### Check the results
str(btc_cleaned)

## Convert data types

### Specify the columns to clean
cols_to_clean <- c("close", "open",
                   "high", "low")

### Remove comma
for (col in cols_to_clean) {
  btc_cleaned[[col]] <- gsub(",",
                             "",
                             btc_cleaned[[col]])
}

### Remove "K"
btc_cleaned$volume <- gsub("K",
                           "",
                           btc_cleaned$volume)

### Remove "%"
btc_cleaned$change <- gsub("%",
                           "",
                           btc_cleaned$change)

### Specify columns to convert to numeric
cols_to_num <- colnames(btc_cleaned[, -1])

### Convert data type to numeric
for (col in cols_to_num) {
  btc_cleaned[[col]] <- as.numeric(btc_cleaned[[col]])
}

### Check the results
str(btc_cleaned)

ผลลัพธ์:

'data.frame':	4955 obs. of  7 variables:
 $ date  : chr  "Feb 09, 2024" "Feb 08, 2024" "Feb 07, 2024" "Feb 06, 2024" ...
 $ close : num  47545 45293 44340 43088 42697 ...
 $ open  : num  45293 44346 43088 42698 42581 ...
 $ high  : num  47710 45579 44368 43376 43532 ...
 $ low   : num  45254 44336 42784 42567 42273 ...
 $ volume: num  86.8 66.4 48.6 33.3 39.3 ...
 $ change: num  4.97 2.15 2.91 0.91 0.27 -0.99 -0.44 0.26 1.18 -0.85 ...

ตอนนี้ ทุก column (ยกเว้น date) ถูกเปลี่ยนให้เป็น numeric เรียบร้อยแล้ว และชุดข้อมูลก็พร้อมจะถูกนำไปใช้งานต่อแล้ว

🔄️ Converting

เรามาดูวิธีแรกในการทำงานกับ time series data กัน นั่นคือ:

การแปลงข้อมูลให้เป็น datetime
การแปลงข้อมูลให้เป็น time series

1️⃣ การแปลงข้อมูลให้เป็น Datetime

ในการทำงานกับ time series data เราจะต้องแปลงข้อมูลประเภทอื่น ๆ เช่น character และ numeric ให้เป็น datetime ก่อน เช่น date (เช่น “Feb 09, 2024”) ใน Bitcoin dataset

เราสามารถแปลงข้อมูลจาก character เป็น Date ได้ 2 วิธี ดังนี้:

วิธีที่ 1. ใช้ as.Date() ซึ่งเป็น base R function และต้องการ input 2 อย่าง ได้แก่:

x: ข้อมูลที่ต้องการแปลงเป็น Date
format: format ของวันเวลาของข้อมูลต้นทาง (เช่น วัน-เดือน-ปี, ปี-เดือน-วัน, เดือน-วัน-ปี)

# Convert `date` to Date
btc_cleaned$date <- as.Date(btc_cleaned$date,
                            format = "%b %d, %Y")

Note: ศึกษาการกำหนด format ได้ที่ strptime: Date-time Conversion Functions to and from Character และ Date Formats in R

วิธีที่ 2. ใช้ as_date() จาก lubridate ซึ่งเป็นที่นิยมมากกว่า as.Date() เพราะ as_date() แปลงข้อมูลให้อัตโนมัติ โดยเราไม่ต้องกำหนด format ให้:

# Convert `date` to Date
btc_cleaned$date <- as_date(btc_cleaned$date)

จากตัวอย่าง จะเห็นได้ว่า as_date() ใช้งานง่ายกว่า as.Date()

ทั้งนี้ ทั้งสองวิธีให้ผลลัพธ์เหมือนกัน คือ แปลง character เป็น Date ซึ่งเราสามารถเช็กผลลัพธ์ได้ด้วย str():

# Check the results
str(btc_cleaned)

ผลลัพธ์:

'data.frame':	4955 obs. of  7 variables:
 $ date  : Date, format: "2024-02-09" "2024-02-08" "2024-02-07" ...
 $ close : num  47545 45293 44340 43088 42697 ...
 $ open  : num  45293 44346 43088 42698 42581 ...
 $ high  : num  47710 45579 44368 43376 43532 ...
 $ low   : num  45254 44336 42784 42567 42273 ...
 $ volume: num  86.8 66.4 48.6 33.3 39.3 ...
 $ change: num  4.97 2.15 2.91 0.91 0.27 -0.99 -0.44 0.26 1.18 -0.85 ...

จะเห็นได้ว่า ตอนนี้ price เปลี่ยนเป็นจาก character เป็น Date เรียบร้อยแล้ว

2️⃣ การแปลงข้อมูลให้เป็น Time Series

นอกจากการแปลงข้อมูลบางส่วนให้เป็น datetime แล้ว เรายังต้องการแปลง dataset จาก data frame ให้เป็น time series เพื่อมี index ของข้อมูล เป็น datetime

ในภาษา R เรามี time series 2 ประเภทหลัก ได้แก่:

ts ของ base R ซึ่งเหมาะกับข้อมูลที่มีความถี่คงที่ (เช่น ข้อมูลมีระยะห่างกัน 1 วัน)
zoo จาก zoo package ซึ่งเหมาะกับข้อมูลที่มีความถี่ไม่คงที่ (เช่น บางข้อมูลห่างกัน 1 วัน และบางข้อมูลห่างกัน 1 เดือน)

ในตัวอย่าง Bitcoin เราจะแปลง dataset ให้เป็น zoo กัน:

# Convert data frame to zoo object

## Convert
btc_zoo <- zoo(btc_cleaned[, -1],
               order.by = btc_cleaned$date)

## Check the results
head(btc_zoo)

ผลลัพธ์:

           close open high low volume change
2010-07-18   0.1  0.0  0.1 0.1   0.08      0
2010-07-19   0.1  0.1  0.1 0.1   0.57      0
2010-07-20   0.1  0.1  0.1 0.1   0.26      0
2010-07-21   0.1  0.1  0.1 0.1   0.58      0
2010-07-22   0.1  0.1  0.1 0.1   2.16      0
2010-07-23   0.1  0.1  0.1 0.1   2.40      0

จะเห็นได้ว่า datetime ถูกเปลี่ยนให้เป็น index ของข้อมูลแล้ว

🤷‍♂️ Missing Values

ในบางครั้ง เราอาจมี time series data ที่มีข้อมูลขาดหายไป ซึ่งเราอาจต้องการแทนที่ด้วยค่าบางอย่าง

ในตัวอย่าง Bitcoin dataset เรามี missing values อยู่ใน volume ซึ่งเราสามารถเช็กได้ด้วย colSums() และ is.na() แบบนี้:

# Count the missing values in each column
colSums(is.na(btc_zoo))

ผลลัพธ์:

 close   open   high    low volume change 
     0      0      0      0    271      0

zoo package มี 3 functions ที่ช่วยให้เราแทนค่า missing values ได้ ดังนี้:

No.	Function	Description
1	`na.fill()`	แทนค่าด้วยค่าที่เรากำหนด
2	`na.locf()`	แทนค่าตามค่าของข้อมูลที่มาก่อนหน้า
3	`na.approx()`	แทนค่าด้วยตามข้อมูลที่มาก่อนและหลัง missing values

สำหรับ Bitcoin dataset เราจะใช้ na.approx() เนื่องจากข้อมูล Bitcoin มีความผันผวนและการแทนค่าแบบยืดหยุ่นจะดีกว่าการแทนค่าแบบกำหนดตายตัว:

# Impute with na.approx()
btc_zoo <- na.approx(btc_zoo)

# Check the results
colSums(is.na(btc_zoo))

ผลลัพธ์:

 close   open   high    low volume change 
     0      0      0      0      0      0

ตอนนี้ missing values ถูกแทนค่าเรียบร้อย

📈 Plotting

การสร้างกราฟสามารถช่วยให้เราเข้าใจ time series data ได้ง่ายขึ้น

เราสามารถทำได้ง่าย ๆ ด้วย autoplot.*() เช่น autoplot.ts() และ autplot.zoo()

ใน Bitcoin dataset เราสามารถสร้างกราฟราคาปิดของ Bitcoin ในช่วง 2010**–**2024 ได้ด้วย autoplot.zoo() แบบนี้:

# Plot the time series
autoplot.zoo(btc_zoo[, "close"]) +
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin Closing Price Over Time (2010–2024)",
       x = "Time",
       y = "Closing Price (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

Note: ในตัวอย่าง เราใช้ ggplot2 ช่วยตกแต่งกราฟให้เข้าใจง่ายขึ้น, สำหรับคนที่สนใจ สามารถศึกษาการใช้ ggplot2 ได้ที่ วิธีใช้ ggplot2 เพื่อสร้างกราฟอย่างมืออาชีพระดับโลก แบบ BBC และ Financial Times ในภาษา R — ตัวอย่างการสำรวจข้อมูลเพนกวินจาก palmerpenguins

✂️ Subsetting

ในการทำงานกับ time series data เราอาจต้องการตัดแบ่งข้อมูลบางส่วนมาเพื่อสำรวจการเปลี่ยนแปลงของข้อมูลในช่วงเวลาที่ต้องการ

ยกตัวอย่างเช่น เราอยากจะสำรวจการเปลี่ยนแปลงราคาปิดของ Bitcoin ในช่วง Crypto Winter (Nov 2021 – Dec 2022)

ในภาษา R เราสามารถตัดแบ่ง time series ได้ 3 วิธี ได้แก่:

window()
Boolean masking
Specific date

1️⃣ window()

window() เป็น base R function ซึ่งช่วยเราเลือกช่วงเวลาได้ด้วย input 3 อย่าง:

x: time series data
start: วันเริ่มต้น (inclusive)
end: วันสิ้นสุด (inclusive)

เราสามารถใช้ window() สร้างกราฟได้แบบนี้:

# Subset
btc_winter_1 <- window(x = btc_zoo,
                       start = as.Date("2021-11-01"),
                       end = as.Date("2022-12-31"))

# Plot the results
autoplot.zoo(btc_winter_1[, "close"]) +
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin Closing Price During the Crypto Winter (2021–2022)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

2️⃣ Boolean Masking

นอกจาก window() แล้ว เราสามารถเลือกช่วงเวลาด้วยการสร้าง Boolean mask ขึ้นมา แบบนี้:

# Create a Boolean masking
crypto_winter_index <- 
  index(btc_zoo) >= "2021-11-01" &
  index(btc_zoo) <= "2022-12-31"

# Subset
btc_winter_2 <- btc_zoo[crypto_winter_index, ]

# Plot the results
autoplot.zoo(btc_winter_2[, "close"]) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin Closing Price During the Crypto Winter (2021–2022)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

จะเห็นว่า กราฟของ Boolean masking เหมือนกับการใช้ window()

3️⃣ Specific Date

ในกรณีที่เราต้องการเลือกวันแทนช่วงเวลา (เช่น ดูข้อมูล Bitcoin เมื่อเกิด halving ครั้งแรก) เราสามารถเขียนภาษา R ได้แบบนี้:

# Subset
btc_first_halving <- btc_zoo["2012-11-28"]

# Print the result
btc_first_halving

ผลลัพธ์:

           close open high  low volume change
2012-11-28  12.4 12.2 12.4 12.1  30.69   1.23

🧮 Aggregating

ในภาษา R เราสามารถสรุปข้อมูล time series data (เช่น หาค่าเฉลี่ย) ได้ด้วย apply.*() จาก zoo package:

Function	Description
`apply.daily()`	สรุปข้อมูลรายวัน
`apply.weekly()`	สรุปข้อมูลรายสัปดาห์
`apply.quarterly()`	สรุปข้อมูลรายไตรมาส
`apply.yearly()`	สรุปข้อมูลรายปี

โดย apply.*() ต้องการ input 2 อย่าง คือ:

x: time series data
FUN: การสรุปข้อมูล (เช่น mean, median, max, min)

ยกตัวอย่างเช่น หาค่าเฉลี่ยราคาปิดของ Bitcoin ในแต่ละปี:

# Example 1. View yearly mean closing price
btc_yr_mean <- apply.yearly(btc_zoo[, "close"],
                            FUN = mean)

# Plot the results
autoplot.zoo(btc_yr_mean) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin Yearly Mean Closing Price (2010–2024)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

หรือหาราคาปิดสูงสุดในแต่ละไตรมาส:

# Example 2. View quarterly max closing price
btc_qtr_max <- apply.quarterly(btc_zoo[, "close"],
                               FUN = max)

# Plot the results
autoplot.zoo(btc_qtr_max) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin Quarterly Maximum Closing Price (2010–2024)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

ทั้งนี้ ถ้าเราต้องการกำหนดความถี่ของเวลาเอง เราสามารถใช้ endpoints() คู่กับ period.apply() ได้

โดย endpoints() ต้องการ input 3 อย่าง:

x: time series data
on: ระดับช่วงเวลา (เช่น ปี, เดือน, สัปดาห์)
k: จำนวน

และ period.apply() ต้องการ input 3 อย่าง คือ:

x: time series data
INDEX: variable ที่เก็บ endpoints() เอาไว้
FUN: function ที่ใช้สรุปข้อมูล (เช่น mean())

ยกตัวอย่างเช่น เราต้องการหาค่าเฉลี่ยราคาปิด Bitcoin ทุก ๆ 3 เดือน เราจะเขียนภาษา R ได้แบบนี้:

# Create 3-month interval
three_month_eps <- endpoints(x = btc_zoo,
                             on = "months",
                             k = 3)

# Apply the interval
btc_three_month_data <- period.apply(x = btc_zoo[, "close"],
                                     INDEX = three_month_eps,
                                     FUN = mean)

# Plot the results
autoplot.zoo(btc_three_month_data) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin 3-Month Average Closing Price (2010–2024)",
       x = "Date",
       y = "Closing Price (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

🪟 Rolling Window

Rolling window (moving window หรือ sliding window) คือ การสรุปข้อมูลตามช่วงเวลาที่กำหนด โดยช่วงเวลาในการคำนวณจะเคลื่อนตัวไปเรื่อย ๆ

เช่น เราต้องการหาค่าเฉลี่ยของราคาปิด Bitcoin ในช่วง 7 วันล่าสุด การใช้ rolling window จะทำให้เราได้ข้อมูลดังนี้:

Date	Average
2024-01-07	ค่าเฉลี่ยระหว่างวันที่ 1 ถึง 7
2024-01-08	ค่าเฉลี่ยระหว่างวันที่ 2 ถึง 8
2024-01-09	ค่าเฉลี่ยระหว่างวันที่ 3 ถึง 9

ในภาษา R เราสามารถสร้าง rolling window ได้ 2 วิธี:

ใช้ roll*()
ใช้ rollapply()

**1️⃣ roll*()**

เราสามารถสร้าง rolling window ได้ด้วย roll*() จาก zoo package เช่น:

Function	Description
`rollmean()`	หาค่าเฉลี่ย
`rollmedian()`	หาค่ากลาง
`rollmax()`	หาค่าสูงสุด
`rollsum()`	หาผลรวม

roll*() ต้องการ input 3 อย่าง ได้แก่:

x: time series data
k: ช่วงเวลา (window)
align: กำหนดว่า ค่าสรุปที่ได้จะอยู่ด้วยซ้าย ขวา หรือกลางช่วงเวลา
fill: จะแทนค่า missing values หรือไม่

ยกตัวอย่างเช่น เราต้องการหาค่าเฉลี่ยราคาปิด Bitcoin ในช่วง 30 วัน เราจะใช้ rollmean() แบบนี้:

# Create the window for mean price
btc_30_days_roll_mean <- rollmean(x = btc_zoo,
                                  k = 30,
                                  align = "right",
                                  fill = NA)

# Plot the results
autoplot.zoo(btc_30_days_roll_mean[, "close"]) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin 30-Day Rolling Mean Price (2010–2024)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

2️⃣ rollapply()

ในกรณีที่เราต้องการกำหนดการสรุปข้อมูลใน rolling window (เช่น หาค่าต่ำสุด ซึ่งไม่มีใน zoo package) เราสามารถทำได้ด้วย rollapply() ซึ่งต้องการ input 4 อย่างดังนี้:

data: time series data
width: ช่วงเวลา (window)
FUN: การสรุปข้อมูล (เช่น min)
align: กำหนดว่า ค่าสรุปที่ได้จะอยู่ด้วยซ้าย ขวา หรือกลางช่วงเวลา
fill: จะแทนค่า missing values หรือไม่

ยกตัวอย่างเช่น หาราคาปิดต่ำสุดของทุก ๆ 30 วัน:

# Creating a rolling window with min() function
btc_30_days_roll_min <- rollapply(data = btc_zoo,
                                  width = 30,
                                  FUN = min,
                                  align = "right",
                                  fill = NA)

# Plot the results
autoplot.zoo(btc_30_days_roll_min[, "close"]) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin 30-Day Rolling Minimum Price (2010–2024)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

➡️ Expanding Window

Expanding window เป็นการสรุปข้อมูลแบบสะสม เช่น:

Date	Average
2024-01-01	หาค่าเฉลี่ยของวันที่ 1
2024-01-02	หาค่าเฉลี่ยของวันที่ 1 + 2
2024-01-03	หาค่าเฉลี่ยของวันที่ 1 + 2 + 3
2024-01-04	หาค่าเฉลี่ยของวันที่ 1 + 2 + 3 + 4
2024-01-05	หาค่าเฉลี่ยของวันที่ 1 + 2 + 3 + 4 + 5

ในภาษา R เราสามารถสร้าง expanding window ได้ด้วย rollapply() โดยเราต้องกำหนดให้:

width เป็นชุดตัวเลข (เช่น 1, 2, 3, …)
align = "right" เสมอ

ยกตัวอย่างเช่น สร้าง expanding window แบบค่าเฉลี่ยราคาปิดของเดือนมกราคม 2024:

# Subset for Jan 2024 data
btc_jan_2024 <- window(x = btc_zoo,
                       start = as.Date("2024-01-01"),
                       end = as.Date("2024-01-31"))

# Create a sequence of widths
btc_jan_2024_width <- seq_along(btc_jan_2024)

# Create an expanding window for mean price
btc_exp_mean <- rollapply(data = btc_2024,
                          width = btc_jan_2024_width,
                          FUN = mean,
                          align = "right",
                          fill = NA)

# Plot the results
autoplot.zoo(btc_exp_mean[, "close"]) +  
  
  ## Adjust line colour
  geom_line(color = "blue") +
  
  ## Add title and labels
  labs(title = "Bitcoin Expanding Mean Price (Jan 2024)",
       x = "Time",
       y = "Closing (USD)") +
  
  ## Use minimal theme
  theme_minimal()

ผลลัพธ์:

😎 Summary

ในบทความนี้ เราได้ทำความรู้จักกับ 7 วิธีการทำงานกับ time series ในภาษา R โดยใช้ base R, lubridate, และ zoo package กัน:

วิธีการทำงาน 1. Converting:

Function	For
`as.Date()`	แปลง character เป็น Date
`as_date()`	แปลง character เป็น Date
`zoo()`	แปลง data frame เป็น `zoo`

วิธีการทำงาน #2. Missing values:

Function	For
`na.fill()`	แทนค่า missing values ด้วยค่าที่กำหนด
`na.locf()`	แทนค่า missing values ตามค่าของข้อมูลที่มาก่อนหน้า
`na.approx()`	แทนค่า missing values ด้วยการวิเคราะห์ข้อมูลที่มาก่อนและหลัง missing values

วิธีการทำงาน #3. Plotting:

Function	For
`autoplot.zoo()`	สร้างกราฟจาก `zoo`

วิธีการทำงาน #4. Subsetting:

Function	For
`window()`	ดึงข้อมูลช่วงเวลาที่ต้องการ
`btc_zoo[boolean_mask, ]`	ดึงข้อมูลช่วงเวลาที่ต้องการ
`btx_zoo[x, y]`	ดึงข้อมูลวันที่ต้องการ

วิธีการทำงาน #5. Aggregating:

Function	For
`apply.daily()`	สรุปข้อมูลรายวัน
`apply.weekly()`	สรุปข้อมูลรายสัปดาห์
`apply.quarterly()`	สรุปข้อมูลรายไตรมาส
`apply.yearly()`	สรุปข้อมูลรายปี
`endpoints()` และ `period.apply()`	สรุปข้อมูลตามช่วงเวลาที่กำหนดเอง

วิธีการทำงาน #6. Rolling window:

Function	For
`rollmean()`	rolling window แบบหาค่าเฉลี่ย
`rollmedian()`	rolling window แบบหาค่ากลาง
`rollmax()`	rolling window แบบหาค่าสูงสุด
`rollsum()`	rolling window แบบหาผลรวม
`rollapply()`	rolling window แบบกำหนดการสรุปข้อมูลเอง

วิธีการทำงาน #7. Expanding window:

Function	For
`seq_along()` และ `rollapply()`	สร้าง expanding window

😺 View Code in GitHub

ดู code และ dataset ในบทความนี้ได้ที่ GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-07-17

สอนปลูกต้นไม้ในภาษา R (ภาค 3): 5 ขั้นตอนการสร้าง Boosted Trees ด้วย xgboost package ในภาษา R — ตัวอย่างการทำนายระดับการกินน้ำมันของรถใน mpg dataset
ในบทความนี้ เราจะไปทำความรู้จักกับ eXtreme Gradient Boosting (XGBoost) และวิธีสร้าง XGBoost model ในภาษา R ด้วย xgboost package กัน

ถ้าพร้อมแล้ว ไปเริ่มกันเลย
อ่านเกี่ยวกับการปลูกต้นไม้ในภาษา R ภาคก่อน ๆ ได้ที่:
- ภาค 1: ปลูกด้วย rpart และ randomForest packages
- ภาค 2: ปลูกด้วย ranger package
🚀 XGBoost คืออะไร?

XGBoost เป็น machine learning model ที่จัดอยู่ในกลุ่ม tree-based models หรือ models ที่ทำนายข้อมูลด้วย decision tree อย่าง single decision tree และ random forest

ใน XGBoost, decision trees จะถูกสร้างขึ้นมาเป็นรอบ ๆ โดยในแต่ละรอบ decision trees ใหม่จะเรียนรู้จากความผิดพลาดของรอบก่อน ซึ่งจะทำให้ decision trees ใหม่มีความสามารถที่ดีขึ้นเรื่อย ๆ

เมื่อสิ้นสุดการสร้าง XGBoost ใช้ผลรวมของ decision trees ทุกต้นในการทำนายข้อมูล ดังนี้:
- Regression problem: หาค่าเฉลี่ยแบบถ่วงน้ำหนักจากทุกต้น
- Classification problem: ทำนายผลลัพธ์ด้วยค่าเฉลี่ยความน่าจะเป็นจากทุกต้น
💻 XGBoost ในภาษา R

ในภาษา R เราสามารถสร้าง XGBoost ได้ด้วย xgboost package ใน 5 ขั้นตอน ได้แก่:
1. Install and load the package
2. Load and prepare the data
3. Split the data
4. Train the model
5. Evaluate the model
1️⃣ Install & Load the Package

ในขั้นแรก ให้เราติดตั้งและเรียกใช้งาน xgboost package

ติดตั้ง:
```
# Install
install.packages("xgboost")
```
เรียกใช้งาน:
```
# Load
library(xgboost)
```
2️⃣ Load & Prepare the Data

ในขั้นตอนที่สอง ให้เราโหลดและเตรียมข้อมูลที่จะใช้สร้าง XGBoost model โดยในบทความนี้ เราจะใช้ mpg dataset จาก ggplot2 package กัน

mpg ประกอบด้วยข้อมูลรถและระดับการใช้น้ำมัน และจุดประสงค์ของเรา คือ ทำนายระดับการกินน้ำมันเมื่อรถวิ่งบน highway (hwy)

เราสามารถโหลด mpg ได้ผ่าน ggplot2:
```
# Install ggplot2
install.packages("ggplot2")

# Load ggplot2
library(ggplot2)

# Load the dataset
data(mpg)
```
เมื่อโหลด dataset แล้ว เราสามารถสำรวจข้อมูลได้ด้วย head():
```
# Preview
head(mpg)
```
ผลลัพธ์:
```
# A tibble: 6 × 11
  manufacturer model displ  year   cyl trans      drv     cty   hwy fl    class  
  <chr>        <chr> <dbl> <int> <int> <chr>      <chr> <int> <int> <chr> <chr>  
1 audi         a4      1.8  1999     4 auto(l5)   f        18    29 p     compact
2 audi         a4      1.8  1999     4 manual(m5) f        21    29 p     compact
3 audi         a4      2    2008     4 manual(m6) f        20    31 p     compact
4 audi         a4      2    2008     4 auto(av)   f        21    30 p     compact
5 audi         a4      2.8  1999     6 auto(l5)   f        16    26 p     compact
6 audi         a4      2.8  1999     6 manual(m5) f        18    26 p     compact
```
และดูโครงสร้างข้อมูลด้วย str():
```
# View the tructure
str(mpg)
```
ผลลัพธ์:
```
tibble [234 × 11] (S3: tbl_df/tbl/data.frame)
 $ manufacturer: chr [1:234] "audi" "audi" "audi" "audi" ...
 $ model       : chr [1:234] "a4" "a4" "a4" "a4" ...
 $ displ       : num [1:234] 1.8 1.8 2 2 2.8 2.8 3.1 1.8 1.8 2 ...
 $ year        : int [1:234] 1999 1999 2008 2008 1999 1999 2008 1999 1999 2008 ...
 $ cyl         : int [1:234] 4 4 4 4 6 6 6 4 4 4 ...
 $ trans       : chr [1:234] "auto(l5)" "manual(m5)" "manual(m6)" "auto(av)" ...
 $ drv         : chr [1:234] "f" "f" "f" "f" ...
 $ cty         : int [1:234] 18 21 20 21 16 18 18 18 16 20 ...
 $ hwy         : int [1:234] 29 29 31 30 26 26 27 26 25 28 ...
 $ fl          : chr [1:234] "p" "p" "p" "p" ...
 $ class       : chr [1:234] "compact" "compact" "compact" "compact" ...
```
จากผลลัพธ์ เราจะเห็นได้ว่า mpg มี columns ที่เราต้องปรับจาก character เป็น factor อยู่ เช่น manufacturer, model ซึ่งเราสามารถปรับได้ดังนี้:
```
# Convert character columns to factor

## Get character columns
chr_cols <- c("manufacturer",
              "model",
              "trans",
              "drv",
              "fl",
              "class")

## For-loop through the character columns
for (col in chr_cols) {
  mpg[[col]] <- as.factor(mpg[[col]])
}

## Check the results
str(mpg)
```
ผลลัพธ์:
```
tibble [234 × 11] (S3: tbl_df/tbl/data.frame)
 $ manufacturer: Factor w/ 15 levels "audi","chevrolet",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ model       : Factor w/ 38 levels "4runner 4wd",..: 2 2 2 2 2 2 2 3 3 3 ...
 $ displ       : num [1:234] 1.8 1.8 2 2 2.8 2.8 3.1 1.8 1.8 2 ...
 $ year        : int [1:234] 1999 1999 2008 2008 1999 1999 2008 1999 1999 2008 ...
 $ cyl         : int [1:234] 4 4 4 4 6 6 6 4 4 4 ...
 $ trans       : Factor w/ 10 levels "auto(av)","auto(l3)",..: 4 9 10 1 4 9 1 9 4 10 ...
 $ drv         : Factor w/ 3 levels "4","f","r": 2 2 2 2 2 2 2 1 1 1 ...
 $ cty         : int [1:234] 18 21 20 21 16 18 18 18 16 20 ...
 $ hwy         : int [1:234] 29 29 31 30 26 26 27 26 25 28 ...
 $ fl          : Factor w/ 5 levels "c","d","e","p",..: 4 4 4 4 4 4 4 4 4 4 ...
 $ class       : Factor w/ 7 levels "2seater","compact",..: 2 2 2 2 2 2 2 2 2 2 ...
```
ตอนนี้ columns ที่เราต้องการถูกเปลี่ยนเป็น factor เรียบร้อยแล้ว

3️⃣ Split the Data

ในขั้นที่สาม เราจะทำ 3 อย่างด้วยกัน คือ:
1. แยกตัวแปรต้น (x) และตัวแปรตาม (y) ออกจากกัน
2. แบ่งข้อมูลออกเป็น training และ test sets
3. แปลงข้อมูลให้เป็น DMatrix
ข้อที่ 1. เราสามารถแยกตัวแปรต้นและตัวแปรตามออกจากกันได้ดังนี้:
```
# Separate the features from the outcome

## Get the features
x <- mpg[, !names(mpg) %in% "hwy"]

## One-hot encode the features
x <- model.matrix(~ . - 1,
                  data = x)

## Get the outcome
y <- mpg$hwy
```
สังเกตว่า ตอนที่เราแยกตัวแปรต้น เราแปลงตัวแปรเหล่านี้เป็น 0, 1 ด้วย one-hot encoding ด้วย เนื่องจาก xgboost ต้องการตัวแปรต้นที่เป็น numeric

.

ข้อที่ 2. จากนั้น เราจะแบ่ง dataset เป็น training (80%) และ test sets (20%) ดังนี้:
```
# Split the data

## Set seed for reproducibility
set.seed(360)

## Get training index
train_index <- sample(1:nrow(x),
                      nrow(x) * 0.8)

## Create x, y train
x_train <- x[train_index, ]
y_train <- y[train_index]

## Create x, y test
x_test <- x[-train_index, ]
y_test <- y[-train_index]

## Check the results
cat("TRAIN SET", "\\n")
cat("1. Data in x_train:", nrow(x_train), "\\n")
cat("2. Data in y_train:", length(y_train), "\\n")
cat("---", "\\n", "TEST SET", "\\n")
cat("1. Data in x_test:", nrow(x_test), "\\n")
cat("2. Data in y_test:", length(y_test), "\\n")
```
ผลลัพธ์:
```
TRAIN SET
1. Data in x_train: 187
2. Data in y_train: 187
---
TEST SET
1. Data in x_test: 47
2. Data in y_test: 47
```
.

ข้อที่ 3. สุดท้าย เราจะแปลง x, y เป็น DMatrix ซึ่งเป็น object ที่ xgboost ใช้ในการสร้าง XGboost model ดังนี้:
```
# Convert to DMatrix

## Training set
train_set <- xgb.DMatrix(data = x_train,
                         label = y_train)

## Test set
test_set <- xgb.DMatrix(data = x_test,
                        label = y_test)

## Check the results
train_set
test_set
```
ผลลัพธ์:
```
TRAIN SET
xgb.DMatrix  dim: 187 x 77  info: label  colnames: yes
---
TEST SET
xgb.DMatrix  dim: 47 x 77  info: label  colnames: yes
```
4️⃣ Train the Model

ในขั้นที่สี่ เราจะสร้าง XGBoost model ด้วย xgb.train() ซึ่งต้องการ 5 arguments ดังนี้:
```
xgb.train(params, data, nrounds, watchlist, verbose)
```
1. params = hyperparametre ที่ต้องการใช้สร้าง model ที่ดีที่สุด
2. data = training set ที่ใช้สร้าง model
3. nrounds = จำนวนครั้งในการในสร้าง decision trees
4. watchlist = ชุดข้อมูลที่ต้องการใช้ประเมินความสามารถของ model
5. verbose = พิมพ์ข้อมูลในระหว่างการสร้าง model (1) หรือไม่ (0)
(Note: ศึกษา argument อื่น ๆ ของ xgb.train() ได้ที่ XGBoost Parameters)

สำหรับบทความนี้ เราจะใช้ xgb.train() ดังนี้:
```
# Train the model

## Set hyperparametres
hp <- list(objective = "reg:squarederror",
           eta = 0.1,
           max_depth = 4,
           eval_metric = c("rmse",
                           "mae"))

## Train
xgb_model <- xgb.train(params = hp,
                       data = train_set,
                       nrounds = 50,
                       watchlist = list(train = train_set,
                                        test = test_set),
                       verbose = 1)
```
ผลลัพธ์:
```
[1]	train-rmse:21.083975	test-rmse:22.739357 
[2]	train-rmse:19.045063	test-rmse:20.598582 
[3]	train-rmse:17.204130	test-rmse:18.713079 
[4]	train-rmse:15.549113	test-rmse:16.974701 
[5]	train-rmse:14.053049	test-rmse:15.453560 
[6]	train-rmse:12.707307	test-rmse:14.097377 
[7]	train-rmse:11.495216	test-rmse:12.877722 
[8]	train-rmse:10.402476	test-rmse:11.767320 
[9]	train-rmse:9.413522	test-rmse:10.740546 
[10]	train-rmse:8.525230	test-rmse:9.863130 
[11]	train-rmse:7.722776	test-rmse:9.068840 
[12]	train-rmse:7.000648	test-rmse:8.357181 
[13]	train-rmse:6.346603	test-rmse:7.687483 
[14]	train-rmse:5.758685	test-rmse:7.091249 
[15]	train-rmse:5.229548	test-rmse:6.557082 
[16]	train-rmse:4.753713	test-rmse:6.079389 
[17]	train-rmse:4.325653	test-rmse:5.651858 
[18]	train-rmse:3.940325	test-rmse:5.275154 
[19]	train-rmse:3.594545	test-rmse:4.938849 
[20]	train-rmse:3.283961	test-rmse:4.627743 
[21]	train-rmse:3.003089	test-rmse:4.352060 
[22]	train-rmse:2.747553	test-rmse:4.110172 
[23]	train-rmse:2.519617	test-rmse:3.889650 
[24]	train-rmse:2.314957	test-rmse:3.691806 
[25]	train-rmse:2.133630	test-rmse:3.499208 
[26]	train-rmse:1.969083	test-rmse:3.330280 
[27]	train-rmse:1.823011	test-rmse:3.181541 
[28]	train-rmse:1.693565	test-rmse:3.045308 
[29]	train-rmse:1.575817	test-rmse:2.919070 
[30]	train-rmse:1.469256	test-rmse:2.812063 
[31]	train-rmse:1.375599	test-rmse:2.700515 
[32]	train-rmse:1.292928	test-rmse:2.615973 
[33]	train-rmse:1.218867	test-rmse:2.541929 
[34]	train-rmse:1.151134	test-rmse:2.462113 
[35]	train-rmse:1.092395	test-rmse:2.404873 
[36]	train-rmse:1.039158	test-rmse:2.336600 
[37]	train-rmse:0.993882	test-rmse:2.291398 
[38]	train-rmse:0.952062	test-rmse:2.236936 
[39]	train-rmse:0.915935	test-rmse:2.198657 
[40]	train-rmse:0.879957	test-rmse:2.152984 
[41]	train-rmse:0.850423	test-rmse:2.102272 
[42]	train-rmse:0.822475	test-rmse:2.054172 
[43]	train-rmse:0.799025	test-rmse:2.011621 
[44]	train-rmse:0.775398	test-rmse:1.971787 
[45]	train-rmse:0.755066	test-rmse:1.933539 
[46]	train-rmse:0.736655	test-rmse:1.900084 
[47]	train-rmse:0.719087	test-rmse:1.870832 
[48]	train-rmse:0.705279	test-rmse:1.853400 
[49]	train-rmse:0.691914	test-rmse:1.834918 
[50]	train-rmse:0.680016	test-rmse:1.825738 
```
จะเห็นได้ว่า model ในแต่ละรอบมี RMSE หรือตัวบ่งชี้ความคลาดเคลื่อน ที่ลดลงเรื่อย ๆ เนื่องจาก model ใหม่เรียนรู้จากความผิดพลาดของ model ก่อนหน้า

หลังจากสร้าง model เสร็จแล้ว เราสามารถดูรายละเอียดของ model ได้แบบนี้:
```
# Print the model
xgb_model
```
ผลลัพธ์:
```
##### xgb.Booster
raw: 62.4 Kb 
call:
  xgb.train(params = hp, data = train_set, nrounds = 50, watchlist = list(train = train_set, 
    test = test_set), verbose = 1)
params (as set within xgb.train):
  objective = "reg:squarederror", eta = "0.1", max_depth = "4", eval_metric = "rmse", validate_parameters = "mae", objective = "TRUE"
xgb.attributes:
  niter
callbacks:
  cb.print.evaluation(period = print_every_n)
  cb.evaluation.log()
# of features: 77 
niter: 50
nfeatures : 77 
evaluation_log:
  iter train_rmse test_rmse
 <num>      <num>     <num>
     1 21.0839746 22.739357
     2 19.0450628 20.598582
   ---        ---       ---
    49  0.6919137  1.834918
    50  0.6800159  1.825738
```
5️⃣ Evaluate the Model

ในขั้นสุดท้าย เราจะประเมินความสามารถของ model ใน 3 ขั้นตอนกัน คือ:
1. ใช้ model ทำนายตัวแปรตาม
2. คำนวณ MAE, RMSE, และ R squared
3. Print MAE, RMSE, และ R squared
.

ข้อที่ 1. ใช้ model ทำนายตัวแปรตาม ด้วย predict():
```
# Make predictions
y_pred <- predict(xgb_model,
                  newdata = x_test)

# Compare predictions to actual outcomes
results <- data.frame(actual = y_test,
                      predicted = y_pred,
                      error = y_test - y_pred)

# Preview the results
head(results, 10)
```
ผลลัพธ์:
```
   actual predicted      error
1      31  27.81219  3.1878090
2      25  25.89449 -0.8944893
3      30  30.13318 -0.1331844
4      29  26.77814  2.2218552
5      24  24.34723 -0.3472347
6      23  23.58175 -0.5817528
7      19  17.81131  1.1886921
8      12  12.32908 -0.3290768
9      12  12.31534 -0.3153391
10     16  16.25793 -0.2579288
```
.

ข้อที่ 2. คำนวณ MAE, RMSE, และ R squared ซึ่งเป็นตัวชี้วัดความสามารถของ regression models:
```
# Calculate MAE
mae <- mean(abs(results$error))

# Calculate RMSE
rmse <- sqrt(mean((results$error)^2))

# Calculate R squared
ss_res <- sum((results$error)^2)
ss_tot <- sum((results$actual - mean(results$actual))^2)
r_squared <- 1 - (ss_res / ss_tot)
```
.

ข้อที่ 3. แสดงผลลัพธ์:
```
# Print the results
cat("MAE:", round(mae, 2), "\\n")
cat("RMSE:", round(rmse, 2), "\\n")
cat("R squared:", round(r_squared, 2), "\\n")
```
ผลลัพธ์:
```
MAE: 1.23
RMSE: 1.83
R squared: 0.93 
```
จะเห็นได้ว่า model ของเราสามารถอธิบายตัวแปรตามได้ถึง 93% (R squared) และมีความคลาดเคลื่อนโดยเฉลี่ย 1.23 miles per gallon (MAE)

💪 Summary

ในบทความนี้ เราได้ไปทำความรู้จักการสร้าง boosted tree ด้วย xgboost package ในภาษา R ซึ่งมีการทำงาน 5 ขั้นตอนกัน:
1. Install and load the package
2. Load and prepare the data
3. Split the data
4. Train the model
5. Evaluate the model
😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References
✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

🙋 ใครที่กำลังเรียนจิตวิทยาหรือทำงานสายจิตวิทยา และเบื่อที่ต้องใช้ software ราคาแพงอย่าง SPSS และ Excel เพื่อทำข้อมูล

💪 ผมขอแนะนำ R Book for Psychologists หนังสือสอนใช้ภาษา R เพื่อการวิเคราะห์ข้อมูลทางจิตวิทยา ที่เขียนมาเพื่อนักจิตวิทยาที่ไม่เคยมีประสบการณ์เขียน code มาก่อน

ในหนังสือ เราจะปูพื้นฐานภาษา R และพาไปดูวิธีวิเคราะห์สถิติที่ใช้บ่อยกัน เช่น:
- Correlation
- t-tests
- ANOVA
- Reliability
- Factor analysis
🚀 เมื่ออ่านและทำตามตัวอย่างใน R Book for Psychologists ทุกคนจะไม่ต้องพึง SPSS และ Excel ในการทำงานอีกต่อไป และสามารถวิเคราะห์ข้อมูลด้วยตัวเองได้ด้วยความมั่นใจ

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists
Share this:
X
Facebook
Like Loading…
2025-07-03

Generalised Linear Model: วิธีใช้ glm() ในภาษา R เพื่อทำนายข้อมูลที่ไม่ปกติ — Linear, Logistic, และ Poisson Regression

ในบทความนี้ เราจะไปทำความรู้จักกับ generalised linear model (GLM) และวิธีทำ GLM ในภาษา R กัน

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

🤔 GLM คืออะไร?

GLM เป็นเทคนิคทางสถิติที่ใช้ทำนายข้อมูลที่มีการกระจายตัวไม่ปกติ (non-normal distribution) เช่น ข้อมูลที่มีผลลัพธ์เพียง 0 และ 1

GLM ทำนายข้อมูลเหล่านี้โดยการต่อยอดจากสมการเส้นตรง (linear model) และมี 3 องค์ประกอบ ได้แก่:

Family: การกระจายตัวของตัวแปรตาม (y)
Linear predictors: สมการเส้นตัวตรงของตัวแปรต้น (x) หรือตัวแปรทำนาย (predictor)
Link function: function ที่เชื่อมตัวแปรต้นกับตัวแปรตามเข้าด้วยกัน

💻 GLM ในภาษา R

ในภาษา R เราสามารถใช้งาน GLM ได้ผ่าน glm() function ซึ่งต้องการข้อมูล 3 อย่าง:

glm(formula, data, family)

formula = ความสัมพันธ์ระหว่างตัวแปรต้นและตัวแปรตาม ในรูปแบบ y ~ x
data = ชุดข้อมูลที่ใช้ในการวิเคราะห์
family = การกระจายตัวของตัวแปรตาม

จะสังเกตว่า glm() ไม่มี parametre สำหรับ link function ทั้งนี้เป็นเพราะ glm() เรียกใช้ link function ให้อัตโนมัติตาม family ที่เรากำหนด

ทั้งนี้ ประเภทข้อมูล, family, และ link function ที่เราสามารถเรียกใช้ glm() ได้มีดังนี้:

Data	family	Link Function
Normal	`gaussian`	`link = "identity”`
Binomial	`binomial`	`link = "logit”`
Poisson	`poisson`	`link = "log”`
Quasi-poisson	`quasipoisson`	`link = "log”`
Gamma	`Gamma`	`link = "inverse”`

เราไปดูตัวอย่างการใช้งาน glm() เพื่อทำนายและแปลผลกัน

☕ ตัวอย่างข้อมูล: coffee_shop

เราจะไปดูตัวอย่างการใช้ glm() เพื่อทำนายข้อมูล 3 ประเภทกัน:

Linear regression
Logistic regression
Poisson regression

โดยเราจะใช้ตัวอย่างเป็นข้อมูลจำลองชื่อ coffee_shop ซึ่งประกอบด้วยข้อมูลการขายจากร้านกาแฟแห่งหนึ่ง และมีรายละเอียดดังนี้:

No.	Column	Description
1	`day`	วันที่
2	`temp`	อุณหภูมิโดยเฉลี่ยของวัน
3	`promo`	เป็นวันที่มีโปรโมชัน (มี, ไม่มี)
4	`weekend`	เป็นวันหยุดสุดสัปดาห์ (วันหยุด, วันธรรมดา)
5	`sales`	ยอดขาย
6	`customers`	จำนวนลูกค้าในแต่ละวัน
7	`sold_out`	ขายหมด (หมด, ไม่หมด)

ก่อนวิเคราะห์ข้อมูล เราจะสร้าง coffee_shop ตามนี้:

# Generate mock coffee shop dataset (15 days)

## Set seed for reproducibility
set.seed(123)

## Generate
coffee_shop <- data.frame(
  
  ## Generate 15 days
  day = 1:15,
  
  ## Generate daily temperature
  temp = round(rnorm(15,
                     mean = 25,
                     sd = 5),
               1),
  
  ## Generate promotion day
  promo = sample(c(0, 1),
                 15,
                 replace = TRUE),
  
  ## Generate weekend
  weekend = c(0, 0, 0, 0, 1, 1, 0, 0, 0, 1, 1, 0, 0, 0, 1),
  
  ## Generate the number of sales
  sales = round(rnorm(15,
                      mean = 300,
                      sd = 50)),
  
  ## Generate the number of daily customers
  customers = rpois(15,
                    lambda = 80),
  
  ## Generate sold-out
  sold_out = sample(c(0, 1),
                    15,
                    replace = TRUE)
)

## Convert binary variables to factors
coffee_shop$promo <- factor(coffee_shop$promo,
                            levels = c(0, 1),
                            labels = c("NoPromo", "Promo"))

coffee_shop$weekend <- factor(coffee_shop$weekend,
                              levels = c(0, 1),
                              labels = c("Weekday", "Weekend"))

coffee_shop$sold_out <- factor(coffee_shop$sold_out,
                               levels = c(0, 1),
                               labels = c("No", "Yes"))

## View the dataset
print(coffee_shop)

ผลลัพธ์:

   day temp   promo weekend sales customers sold_out
1    1 22.2 NoPromo Weekday   246        73      Yes
2    2 23.8   Promo Weekday   296        76      Yes
3    3 32.8 NoPromo Weekday   354        83       No
4    4 25.4   Promo Weekday   293        87       No
5    5 25.6   Promo Weekend   242        78      Yes
6    6 33.6 NoPromo Weekend   259        91      Yes
7    7 27.3 NoPromo Weekday   334        71      Yes
8    8 18.7 NoPromo Weekday   284        76       No
9    9 21.6 NoPromo Weekday   234        72      Yes
10  10 22.8   Promo Weekend   270        79       No
11  11 31.1 NoPromo Weekend   294        82      Yes
12  12 26.8   Promo Weekday   344        79       No
13  13 27.0   Promo Weekday   292        79       No
14  14 25.6 NoPromo Weekday   316        92       No
15  15 22.2 NoPromo Weekend   139        77      Yes

เราไปดูวิธีทำนายข้อมูลกัน

1️⃣ Linear Regression

Linear regression เป็นการทำนายข้อมูล numeric เช่น ยอดขาย (sales) ซึ่งเราสามารถใช้ glm() ทำนายได้ดังนี้:

# Create a regression model with glm()
linear_reg <- glm(sales ~ temp + promo + weekend,
                  data = coffee_shop,
                  family = gaussian)

เราสามารถดู model ได้ด้วย summary():

# Get model summary
summary(linear_reg)

ผลลัพธ์:

Call:
glm(formula = sales ~ temp + promo + weekend, family = gaussian, 
    data = coffee_shop)

Coefficients:
               Estimate Std. Error t value Pr(>|t|)   
(Intercept)      96.599     59.669   1.619  0.13375   
temp              7.703      2.283   3.373  0.00621 **
promoPromo       23.014     18.537   1.241  0.24025   
weekendWeekend  -73.444     19.654  -3.737  0.00328 **
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for gaussian family taken to be 1222.243)

    Null deviance: 39702  on 14  degrees of freedom
Residual deviance: 13445  on 11  degrees of freedom
AIC: 154.54

Number of Fisher Scoring iterations: 2

จากผลลัพธ์ เราจะเห็นความสำคัญของตัวแปรต้นและ coefficient ซึ่งระบุการเปลี่ยนแปลงของตัวแปรตามการเปลี่ยนของตัวแปรต้น:

ตัวแปรที่สามารถทำนาย sales ได้อย่างมีนัยสำคัญ คือ temp และ weekend (สังเกตจาก **)
promo ไม่สามารถทำนาย sales ได้อย่างมีนัยสำคัญ
Coefficient ของ temp คือ 7.70 ซึ่งหมายถึง อุณหภูมิเปลี่ยน 1 หน่วย ยอดขายจะเพิ่มขึ้น 7.70 หน่วย
Coefficient ของ weekend คือ -73.44 ซึ่งหมายถึง วันหยุดสุดสัปดาห์ ยอดขายจะลดลง 73.44 หน่วย

2️⃣ Logistic Regression

Logistic regression เป็นการทำนายข้อมูลที่เป็นมีผลลัพธ์เพียง 2 ค่า เช่น:

ใช่, ไม่ใช่
ผ่าน, ไม่ผ่าน
ตรง, ไม่ตรง

ใน coffee_shop เรามี sold_out ซึ่งเราสามารถทำนายด้วย logistic regression ด้วย glm() ได้แบบนี้:

# Create a logistic regression model
log_reg <- glm(sold_out ~ temp + promo + weekend,
               data = coffee_shop,
               family = binomial)

จากนั้น ดู model ด้วย summary():

# Get model summary
summary(log_reg)

ผลลัพธ์:

Call:
glm(formula = sold_out ~ temp + promo + weekend, family = binomial, 
    data = coffee_shop)

Coefficients:
               Estimate Std. Error z value Pr(>|z|)
(Intercept)     1.30934    3.92442   0.334    0.739
temp           -0.04448    0.15189  -0.293    0.770
promoPromo     -1.73681    1.31881  -1.317    0.188
weekendWeekend  2.15038    1.45643   1.476    0.140

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 20.728  on 14  degrees of freedom
Residual deviance: 16.426  on 11  degrees of freedom
AIC: 24.426

Number of Fisher Scoring iterations: 4

จะเห็นได้ว่า หน้าตาผลลัพธ์คล้ายกับ linear regression แต่สิ่งที่แตกต่างกัน คือ coefficient อยู่ในรูป log-odd ซึ่งเราสามารถถอดรูปได้ด้วย exp() เพื่อแปลผล:

# Transform coefficient
exp(coef(log_reg))

ผลลัพธ์:

   (Intercept)           temp     promoPromo weekendWeekend 
     3.7037124      0.9564938      0.1760812      8.5881237

เราสามารถแปลผลได้ดังนี้:

Predictor	Coefficient	Interpretation
`temp`	0.96	เมื่ออุณหภูมิสูงขึ้น 1 หน่วย ร้านมีโอกาสขายหมดเพิ่มขึ้น 0.96
`promo`	0.18	เมื่อมีโปรโมชัน ร้านมีโอกาสขายหมดเพิ่มขึ้น 0.18
`weekend`	8.59	เมื่อเป็นวันสุดสัปดาห์ ร้านมีโอกาสขายหมดเพิ่มขึ้น 8.59

3️⃣ Poisson Regression

Poisson regression เป็นการทำนายข้อมูลการนับ (count data) หรือข้อมูลที่เกิดขึ้นในช่วงเวลาที่กำหนด เช่น:

จำนวนรถบนถนนในแต่ละชั่วโมง
จำนวนข้อความที่ได้รับใน 1 วัน
จำนวนสินค้าที่ขายได้ใน 3 เดือน

ใน coffee_shop เรามี customers ซึ่งสามารถทำนาย poisson regression ผ่าน glm() ได้ดังนี้:

# Create a poisson regression model
poisson_reg <- glm(customers ~ temp + promo + weekend,
                   data = coffee_shop,
                   family = poisson)

ดู model:

# Get model summary
summary(poisson_reg)

ผลลัพธ์:

Call:
glm(formula = customers ~ temp + promo + weekend, family = poisson, 
    data = coffee_shop)

Coefficients:
               Estimate Std. Error z value Pr(>|z|)    
(Intercept)    4.108939   0.192499  21.345   <2e-16 ***
temp           0.010077   0.007299   1.381    0.167    
promoPromo     0.010327   0.059616   0.173    0.862    
weekendWeekend 0.012692   0.062739   0.202    0.840    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

(Dispersion parameter for poisson family taken to be 1)

    Null deviance: 7.0199  on 14  degrees of freedom
Residual deviance: 4.8335  on 11  degrees of freedom
AIC: 106.06

Number of Fisher Scoring iterations: 3

Coefficient ของ poisson regression อยู่ในรูป log เช่นเดียวกับ logistic regression ดังนั้น เราจะถอดรูปด้วย exp() ก่อนแปลผล:

# Transform the coefficients
exp(coef(poisson_reg))

ผลลัพธ์:

   (Intercept)           temp     promoPromo weekendWeekend 
     60.882073       1.010127       1.010381       1.012773

จะเห็นได้ว่า coefficient ของทั้ง 3 ตัวแปรต้นอยู่ที่ 1.01 ซึ่งหมายถึง การเปลี่ยนแปลงตัวแปรต้นตัวใดตัวหนึ่ง ทำให้จำนวนลูกค้าเพิ่มขึ้น 1 คน

💪 Summary

ในบทความนี้ เราได้ไปทำความรู้จักกับ GLM ซึ่งเป็นเทคนิคทางสถิติที่ใช้ทำนายข้อมูลที่ไม่ปกติ และได้ดูวิธีการใช้ glm() function ในภาษา R เพื่อทำนายข้อมูล 3 ประเภท:

Linear regression
Logistic regression
Poisson regression

😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References

What is GLM?

GLM in R:

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-06-26

tidymodels: แนะนำวิธีใช้แพ็กเกจ machine learning ที่ทรงพลัง ครบจบ และทันสมัยในภาษา R เพื่อสร้าง ประเมิน และปรับทูน model — ตัวอย่างการทำนายราคาบ้านใน Boston dataset

ในการทำ machine learning (ML) ในภาษา R เรามี packages และ functions ที่หลากหลายให้เลือกใช้งาน ซึ่งแต่ละ package และ function มีวิธีใช้งานที่แตกต่างกันไป

ยกตัวอย่างเช่น:

glm() จาก base R สำหรับสร้าง regression models ต้องการ input 3 อย่าง คือ formula, data, และ family:

glm(formula, data, family)

knn() จาก class package สำหรับสร้าง KNN model ต้องการ input 4 อย่าง คือ ตัวแปรต้นของ training set, ตัวแปรต้นของ test set, ตัวแปรตามของ training set, และค่า k:

knn(train_x, test_x, train_y, k)

rpart() จาก rpart package สำหรับสร้าง decision tree model ต้องการ input 2 อย่าง คือ formula และ data:

rpart(formula, data)

…

การใช้งาน function ที่แตกต่างกันทำให้การสร้าง ML models เกิดความซับซ้อนโดยไม่จำเป็น และทำให้เกิดความผิดพลาดในการทำงานได้ง่าย

tidymodels เป็น package ที่ถูกออกแบบมาเพื่อแก้ปัญหานี้โดยเฉพาะ

tidymodels เป็น meta-package หรือ package ที่รวบรวม packages อื่นเอาไว้ เมื่อเราโหลด tidymodels เราจะสามารถใช้งาน 8 packages ที่ออกแบบมาให้ทำงานร่วมกัน ช่วยให้เราทำงาน ML ได้ครบ loop

ทั้ง 8 packages ใน tidymodels ได้แก่:

No.	Package	ML Phase	For
1	`rsample`	Pre-processing	Data resampling
2	`recipes`	Pre-processing	Feature engineering
3	`parsnip`	Modelling	Model fitting
4	`tune`	Post-processing	Hyperparameter tuning
5	`dials`	Post-processing	Hyperparameter tuning
6	`yardstick`	Post-processing	Model evaluation
7	`broom`	Post-processing	Format ผลลัพธ์ให้ดูง่าย
8	`workflow`	All	รวม pre-processing, modeling, and post-processing ให้เป็น pipeline เดียวกัน

ในบทความนี้ เราจะมาดูวิธีใช้ tidymodels เพื่อสร้าง ML models กัน

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

🔢 Dataset

🏠 Boston

ในบทความนี้ เราจะใช้ Boston dataset จาก MASS package เป็นตัวอย่างในการทำงานกับ tidymodels กัน

Boston เป็นชุดข้อมูลบ้านในเมืองบอสตัน รัฐแมสซาชูเซส ประเทศอเมริกา และมีข้อมูลทั้งหมด 14 columns ดังนี้:

No.	Column	Description
1	`crim`	ระดับอาชญากรรมในแต่ละเขต
2	`zn`	สัดส่วนพื้นที่อาศัย
3	`indus`	สัดส่วนธุรกิจที่เป็น non-retail ในแต่ละเขต
4	`chas`	เป็นพื้นที่ติดกับ Charles River ไหม (1 = ติด, 0 = ไม่ติด)
5	`nox`	ระดับ nitrogen oxide
6	`rm`	จำนวนห้องโดยเฉลี่ย
7	`age`	สัดส่วย unit ที่มีคนเข้าอยู่ ซึ่งถูกสร้างก่อนปี ค.ศ. 1940
8	`dis`	ระยะทางจากพื้นที่ทำงานในเมืองบอสตัน
9	`rad`	ระดับการเข้าถึง radial highways
10	`tax`	ภาษีโรงเรือน
11	`ptratio`	สัดส่วนนักเรียนต่อครูในแต่ละเขต
12	`black`	สัดส่วนผู้อยู่อาศัยที่เป็นคนผิวดำ
13	`lstat`	สัดส่วนของประชากรที่มีฐานะยากจน
14	`medv`	ราคากลางของบ้านที่มีผู้อยู่อาศัย

Note: อ่านรายละเอียดเพิ่มเติมเกี่ยวกับ Boston dataset ได้ที่ A Complete Guide to the Boston Dataset in R

จุดประสงค์ของเราในการทำงานกับ Boston dataset คือ ทำนายราคาบ้าน (medv)

⬇️ Get Boston

ในการใช้งาน Boston เราสามารถเรียกใช้งานได้ผ่าน MASS package ดังนี้:

ติดตั้ง MASS package
โหลด MASS package
โหลด Boston

# Install package
install.packages("MASS")

# Load package
library(MASS)

# Load dataset
data(Boston)

🧹 Prepare Boston

ก่อนเริ่มใช้งาน Boston เราจะต้องปรับประเภทข้อมูลของ chas จาก numeric เป็น factor ก่อน ดังนี้:

# Create a new dataset
bt <- Boston

# Convert `chas` to factor
bt$chas <- factor(bt$chas,
                  levels = c(1, 0),
                  labels = c("tract bounds river", "otherwise"))

# Preview
head(Boston)

ผลลัพธ์:

     crim zn indus chas   nox    rm  age    dis rad tax ptratio  black lstat medv
1 0.00632 18  2.31    0 0.538 6.575 65.2 4.0900   1 296    15.3 396.90  4.98 24.0
2 0.02731  0  7.07    0 0.469 6.421 78.9 4.9671   2 242    17.8 396.90  9.14 21.6
3 0.02729  0  7.07    0 0.469 7.185 61.1 4.9671   2 242    17.8 392.83  4.03 34.7
4 0.03237  0  2.18    0 0.458 6.998 45.8 6.0622   3 222    18.7 394.63  2.94 33.4
5 0.06905  0  2.18    0 0.458 7.147 54.2 6.0622   3 222    18.7 396.90  5.33 36.2
6 0.02985  0  2.18    0 0.458 6.430 58.7 6.0622   3 222    18.7 394.12  5.21 28.7

เราสามารถเช็กผลลัพธ์ได้ด้วย str():

# Check results
str(bt)

ผลลัพธ์:

'data.frame':	506 obs. of  14 variables:
 $ crim   : num  0.00632 0.02731 0.02729 0.03237 0.06905 ...
 $ zn     : num  18 0 0 0 0 0 12.5 12.5 12.5 12.5 ...
 $ indus  : num  2.31 7.07 7.07 2.18 2.18 2.18 7.87 7.87 7.87 7.87 ...
 $ chas   : Factor w/ 2 levels "tract bounds river",..: 2 2 2 2 2 2 2 2 2 2 ...
 $ nox    : num  0.538 0.469 0.469 0.458 0.458 0.458 0.524 0.524 0.524 0.524 ...
 $ rm     : num  6.58 6.42 7.18 7 7.15 ...
 $ age    : num  65.2 78.9 61.1 45.8 54.2 58.7 66.6 96.1 100 85.9 ...
 $ dis    : num  4.09 4.97 4.97 6.06 6.06 ...
 $ rad    : int  1 2 2 3 3 3 5 5 5 5 ...
 $ tax    : num  296 242 242 222 222 222 311 311 311 311 ...
 $ ptratio: num  15.3 17.8 17.8 18.7 18.7 18.7 15.2 15.2 15.2 15.2 ...
 $ black  : num  397 397 393 395 397 ...
 $ lstat  : num  4.98 9.14 4.03 2.94 5.33 ...
 $ medv   : num  24 21.6 34.7 33.4 36.2 28.7 22.9 27.1 16.5 18.9 ...

ตอนนี้ ทุก column มีข้อมูลถูกประเภทแล้ว เราพร้อมที่จะใช้ Boston กับ tidymodels แล้ว

🛩️ tidymodels

🏁 Getting Started

ในการเริ่มใช้งาน tidymodels เราต้องติดตั้งและโหลด package ก่อน

ติดตั้ง (ทำครั้งแรกครั้งเดียว):

# Install
install.packages("tidymodels")

โหลด (ทำทุกครั้งที่เริ่ม session ใหม่):

# Load
library(tidymodels)

🌊 Flows

การใช้งาน tidymodels แบ่งเป็น 2 แบบ:

Standard flow: เขียน code ยาวกว่า เหมาะกับการทำงานที่เราต้องการควบคุมการทำงานแต่ละขั้นด้วยตัวเอง
Workflow: เขียน code สั้นกว่า เหมาะกับการสร้าง ML model อย่างรวดเร็ว

เราไปดูวิธีเขียนแต่ละแบบกัน

🥐 Method #1. Standard Flow

สำหรับการใช้งาน tidymodels แบบ standard flow มีทั้งหมด 7 ขั้นตอน ได้แก่:

Split the data
Create a recipe
Prep and bake
Instantiate a model
Fit the model
Make predictions
Evaluate the model

1️⃣ Split the Data

ในขั้นแรก เราจะแบ่ง dataset ออกเป็น 2 ชุด ได้แก่:

Training set สำหรับสร้าง model
Test set สำหรับทดสอบ model

โดยเราจะใช้ 3 functions จาก tidymodels ช่วย ได้แก่:

No.	Function	For
1	`initial_split()`	กำหนดการแบ่ง dataset
2	`training()`	สร้าง training set
3	`testing()`	สร้าง test set

ซึ่งเราเรียกใช้งานได้ดังนี้:

# Set seed for reproducibility
set.seed(2025)

# Define the training set index
bt_split <- initial_split(data = bt,
                          prop = 0.8,
                          strata = medv)

# Create the training set
bt_train <- training(bt_split)

# Create the test set
bt_test <- testing(bt_split)

2️⃣ Create a Recipe

ในขั้นที่สอง เราจะสร้าง recipe หรือสูตรในการเตรียมข้อมูลกัน

ในขั้นนี้ tidymodels มี 2 functions หลักให้เราใช้งาน ได้แก่:

No.	Function	For
1	`recipe()`	กำหนดความสัมพันธ์ระหว่างตัวแปรต้นและตัวแปรตาม
2	`step_*()`	data cleaning และ feature engineering

ซึ่งเราเรียกใช้งานได้ดังนี้:

# Create a recipe
rec <- recipe(medv ~ .,
              data = bt_train) |>
  
  ## Remove near-zero variance predictors
  step_nzv(all_numeric_predictors()) |>
  
  ## Handle multicollinearity
  step_corr(all_numeric_predictors(),
            threshold = 0.8)

ในตัวอย่าง เราใช้:

recipe() กำหนดตัวแปรตาม (medv) และตัวแปรต้น (ตัวแปรที่เหลือทั้งหมด)
step_nzv() ลบตัวแปรต้นที่มี variance เข้าใกล้ 0
step_corr() จัดการกับ multicollinearity

ศึกษาการใช้งาน step_*() อื่น ๆ ในการเตรียมข้อมูลได้ที่ tidymodels.org – Search recipe steps

3️⃣ Prep & Bake

ในขั้นที่ 3 เราจะจัดเตรียมข้อมูลตาม recipe ที่กำหนด โดยเราจะใช้งาน 2 functions ได้แก่:

No.	Function	For
1	`prep()`	เตรียม recipe
2	`bake()`	เตรียมข้อมูลตาม recipe

โดยเราแบ่งการทำงานเป็น 3 ขั้นตอน คือ:

เตรียม recipe ด้วย prep()
เตรียม training data ด้วย bake()
เตรียม test data ด้วย bake()

# Prepare the recipe
rec_prep <- prep(rec,
                 data = bt_train)

# Bake the training set
bt_train_baked <- bake(rec_prep,
                       new_data = NULL)

# Bake the test set
bt_test_baked <- bake(rec_prep,
                      new_data = bt_test)

4️⃣ Instantiate a Model

ในขั้นที่ 4 เราจะเรียกใช้ algorithm สำหรับ model ของเรา โดยในตัวอย่าง เราจะลองสร้าง decision tree กัน

ในขั้นนี้ เรามี 3 functions จะเรียกใช้งาน ได้แก่:

No.	Function	For
1	`decision_tree()`	สร้าง decision tree *
2	`set_engine()`	กำหนด engine หรือ package ที่ใช้สร้าง model
3	`set_mode()`	กำหนดประเภท model (classification หรือ regression)

Function นี้จะเปลี่ยนตาม model ที่ต้องการ โดยเราสามารถค้นหา model ที่ต้องการได้ที่ tidymodels.org – Search parsnip models

# Instantiate the model
dt_mod <- decision_tree() |>
  
  # Set the engine
  set_engine("rpart") |>
  
  # Set the mode
  set_mode("regression")

ในตัวอย่าง เราใช้:

set_engine() เลือก rpart เป็น engine ในการสร้าง decision tree
set_mode() กำหนด mode เป็น regression เพราะเราต้องการทำนายราคาบ้านซึ่งเป็น continuous variable

5️⃣ Fit the Model

ในขั้นที่ 5 เราจะ train model ด้วย training set ผ่าน fit():

dt_mod_fit <- fit(dt_mod,
                  medv ~ .,
                  data = bt_train_baked)

ตอนนี้ เราก็จะได้ model ที่พร้อมใช้งานมาแล้ว

6️⃣ Make Predictions

ในขั้นที่ 6 เราจะใช้ model ทำนายข้อมูลเพื่อนำไปทดสอบความสามารถในขั้นที่ 7 ต่อไป

ในขั้นนี้ เราจะใช้ predict() ร่วมกับ bind_cols() เพื่อเก็บผลลัพธ์การทำนายเอาไว้:

# Make predictions
dt_results <- predict(dt_mod_fit,
                      new_data = bt_test_baked,
                      type = "numeric") |>
  bind_cols(actual = bt_test_baked$medv)

# Print the results
dt_results

ผลลัพธ์:

# A tibble: 103 × 2
   .pred actual
   <dbl>  <dbl>
 1  34.6   34.7
 2  34.6   33.4
 3  11.8   16.5
 4  17.2   15.6
 5  27.3   30.8
 6  21.7   25  
 7  34.6   35.4
 8  21.7   21.2
 9  27.3   23.9
10  44.3   43.8
# ℹ 93 more rows
# ℹ Use `print(n = ...)` to see more rows

7️⃣ Evaluate the Model

ในขั้นสุดท้าย เราจะวิเคราะห์ความสามารถของ model กัน

tidymodels มี functions สำหรับคำนวณค่าตัวชี้วัดต่าง ๆ เช่น:

Function	For
`accuracy()`	ความแม่นยำในการทำนาย
`roc_auc()`	ความสมดุลในการทำนาย

Note: ศึกษา functions ทั้งหมดได้ที่ tidymodels.org – Metric types

สำหรับบทความนี้ เราจะเลือกใช้ 2 functions ได้แก่:

No.	Function	For
1	`mae()`	Mean absolute error
2	`rmse()`	Root mean squared error

โดยเราเรียกใช้งาน functions ได้ 2 แบบ ดังนี้:

แบบที่ 1. เรียกใช้ด้วยตัวเอง:

# Calculate MAE
dt_mae <- mae(dt_results,
              truth = actual,
              estimate = .pred)

# Calculate RMSE
dt_rmse <- rmse(dt_results,
                truth = actual,
                estimate = .pred)

# Print MAE and RMSE
cat("MAE:", round(dt_mae$.estimate, 2), "\\n")
cat("RMSE:", round(dt_rmse$.estimate, 2), "\\n")

ผลลัพธ์:

MAE:
# A tibble: 1 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 mae     standard        3.70
------------------------------------------ 
RMSE:
# A tibble: 1 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 rmse    standard        5.13

แบบที่ 2. เรียกใช้ผ่าน metric_set() ที่จะรวม functions ไว้ด้วยกัน:

# Define a custom metrics
dt_metrics <- metric_set(mae,
                         rmse)

# Evaluate the model
dt_eva_results <- dt_metrics(dt_results,
                             truth = actual,
                             estimate = .pred)

# Print the results
dt_eva_results

ผลลัพธ์:

# A tibble: 2 × 3
  .metric .estimator .estimate
  <chr>   <chr>          <dbl>
1 mae     standard        3.70
2 rmse    standard        5.13

จะเห็นได้ว่า แบบที่ 2 สะดวกกว่า เพราะเราสามารถคำนวณค่าตัวชี้วัดหลายตัวได้ในครั้งเดียวกัน

🍰 Method #2. Workflow

เราได้ดูวิธีใช้งานแบบ standard flow ไปแล้ว เรามาดูวิธีใช้ tidymodels แบบ workflow ซึ่งมี 6 ขั้นตอน ดังนี้กัน:

Split the data
Create a recipe
Instantiate the model
Bundle the recipe and the model
Fit the model
Evaluate the model

1️⃣ Split the Data

ในขั้นแรก เราจะแบ่ง dataset ออกเป็น 2 ชุด (เหมือนกับ standard flow):

# Set seed for reproducibility
set.seed(2025)

# Define the training set index
bt_split <- initial_split(data = bt,
                          prop = 0.8,
                          strata = medv)

# Create the training set
bt_train <- training(bt_split)

จะสังเกตว่า เราจะไม่ได้สร้าง test set ในครั้งนี้

2️⃣ Create a Recipe

ในขั้นที่ 2 เราจะสร้าง recipe (เหมือนกับ standard flow):

# Create a recipe
rec <- recipe(medv ~ .,
              data = bt_train) |>
  
  ## Remove near-zero variance predictors
  step_nzv(all_numeric_predictors()) |>
  
  ## Handle multicollinearity
  step_corr(all_numeric_predictors(),
            threshold = 0.8)

3️⃣ Instantiate the Model

ในขั้นที่ 3 ของ standard flow เราจะเตรียม recipe และข้อมูลกัน

แต่ใน workflow เราจะสร้าง model (ซึ่งเป็นขั้นที่ 4 ของ standard flow) แทน:

# Instantiate the model
dt_mod <- decision_tree() |>
  
  ## Set the engine
  set_engine("rpart") |>
  
  ## Set the mode
  set_mode("regression")

4️⃣ Bundle the Recipe and the Model

ในขั้นที่ 4 เราจะรวม recipe และ model เข้าด้วยกัน เพื่อทำให้การทำงานต่อจากนี้ง่ายขึ้น ผ่านการใช้ 3 functions ดังนี้:

No.	Function	For
1	`workflow()`	สร้าง `workflow` object
2	`add_recipe()`	เพิ่ม recipe ใน `workflow` object
3	`add_model()`	เพิ่ม model ใน `workflow` object

Note: อ่านเพิ่มเติมเกี่ยวกับ workflow object ได้ที่ tidymodels.org – workflows และ Tidy Modeling with R – A Model Workflow

# Bundle the recipe and the model
dt_wfl <- workflow() |>
  
  ## Add recipe
  add_recipe(rec) |>
  
  ## Add model
  add_model(dt_mod)

5️⃣ Fit the Model

ในขั้นที่ 5 เราจะ train model ด้วย training set โดยเราจะใช้ last_fit() แทน fit()

ทั้งนี้ last_fit() ต้องการ input 3 อย่าง ได้แก่:

No.	Input	Description
1	object	`workflow` object
2	split	object ที่ได้จาก `initial_split()`
3	metrics	ตัวชี้วัดที่เก็บไว้ใน `metric_set()`

# Fit the model
dt_last_fit <- last_fit(dt_wfl,
                        split = bt_split,
                        metrics = metric_set(mae, rmse))

6️⃣ Evaluate the Model

นอกจาก train model แล้ว last_fit() ยังทำหน้าที่อีก 2 อย่าง ได้แก่:

ทำนายข้อมูลจาก model ที่ได้
ประเมิน model ตามตัวชี้วัดที่กำหนดใน metric_set()

นั่นหมายคงามว่า last_fit() ได้ให้ผลลัพธ์ในการประเมิน model มาแล้ว เราเพียงแค่ต้องเรียกผลลัพธ์ออกมาแสดงเท่านั้น ซึ่งเราสามารถทำได้ด้วย 2 functions นี้:

No.	Function	For
1	`collect_predictions()`	ดึงข้อมูลที่ model ทำนาย
2	`collect_metrics()`	ดึงค่าตัวชี้วัดความสามารถของ model

ดึงผลลัพธ์ในการทำนายได้ด้วย collect_predictions():

# Collect predictions
dt_predictions <- collect_predictions(dt_last_fit)

# Print predictions
dt_predictions

ผลลัพธ์:

# A tibble: 103 × 5
   .pred id                .row  medv .config             
   <dbl> <chr>            <int> <dbl> <chr>               
 1  34.6 train/test split     3  34.7 Preprocessor1_Model1
 2  34.6 train/test split     4  33.4 Preprocessor1_Model1
 3  11.8 train/test split     9  16.5 Preprocessor1_Model1
 4  17.2 train/test split    25  15.6 Preprocessor1_Model1
 5  27.3 train/test split    40  30.8 Preprocessor1_Model1
 6  21.7 train/test split    53  25   Preprocessor1_Model1
 7  34.6 train/test split    56  35.4 Preprocessor1_Model1
 8  21.7 train/test split    79  21.2 Preprocessor1_Model1
 9  27.3 train/test split    82  23.9 Preprocessor1_Model1
10  44.3 train/test split    99  43.8 Preprocessor1_Model1
# ℹ 93 more rows
# ℹ Use `print(n = ...)` to see more rows

Note:

.pred คือ ราคาที่ทำนาย
medv คือ ราคาจริง

ดึงตัวบ่งชี้ที่ได้จาก last_fit() ด้วย collect_metrics():

# Collect metrics
dt_metrics <- collect_metrics(dt_last_fit)

# Print metrics
dt_metrics

ผลลัพธ์:

# A tibble: 2 × 4
  .metric .estimator .estimate .config             
  <chr>   <chr>          <dbl> <chr>               
1 mae     standard        3.70 Preprocessor1_Model1
2 rmse    standard        5.13 Preprocessor1_Model1

จะเห็นได้ว่า แม้ workflow จะมีจำนวนขั้นตอนใกล้เคียงกับ standard flow แต่การทำงานแบบ workflow มีความสะดวกกว่ามาก

🍩 Bonus: Hyperparametre Tuning

ส่งท้าย เรามาดูวิธีใช้ tidymodels เพื่อทำ hyperparametre tuning ใน 5 ขั้นตอนกัน:

Prepare
Tune
Select
Fit
Evaluate

1️⃣ Prepare

ในการทำ hyperparametre tuning เราสามารถเขียนได้ทั้งแบบ standard flow และ workflow

ในบทความนี้ เราจะดูวิธีทำแบบ workflow ซึ่งเป็นวิธีที่ง่ายกว่ากัน

ในขั้นแรก เราต้องเตรียมของสำหรับ hyperparametre tuning 4 ข้อ ได้แก่:

No.	Input	Description
1	object	`workflow` object ที่รวม recipe และ model เอาไว้
2	resamples	การแบ่ง training set เพื่อทดสอบ hyperparametres
3	grid	การจับคู่ hyperparametres
4	metrics	ตัวชี้วัดที่ต้องการใช้ประเมิน hyperparametres

ข้อที่ 1. เตรียม object

เราเริ่มจากแบ่งข้อมูล:

# Set seed for reproducibility
set.seed(2025)

# Define the training set index
bt_split <- initial_split(data = bt,
                          prop = 0.8,
                          strata = medv)

# Create the training set
bt_train <- training(bt_split)

สร้าง recipe:

# Create a recipe
rec <- recipe(medv ~ .,
              data = bt_train) |>
  
  ## Remove near-zero variance predictors
  step_nzv(all_numeric_predictors()) |>
  
  ## Handle multicollinearity
  step_corr(all_numeric_predictors(),
            threshold = 0.8)

จากนั้น เรียกใช้ model โดยเราจะต้องกำหนด hyperparametre ที่ต้องการปรับ ด้วย tune():

# Define the tuning parameters
dt_model_tune <- decision_tree(cost_complexity = tune(),
                               tree_depth = tune(),
                               min_n = tune()) |>
  
  ## Set engine
  set_engine("rpart") |>
  
  ## Set mode
  set_mode("regression")

สุดท้าย เรารวม recipe และ model ไว้ใน workflow object:

# Define the workflow with tuning
bt_wfl_tune <- workflow() |>
  
  ## Add recipe
  add_recipe(rec) |>
  
  ## Add model
  add_model(dt_model_tune)

ข้อที่ 2. Resamples ซึ่งในที่นี้ เราจะใช้ k-fold cross-validation ที่แบ่ง training set ออกเป็น 5 ส่วน (4 ส่วนเพื่อ train และ 1 ส่วนเพื่อ test) แบบนี้:

# Set k-fold cross-validation for tuning
hpt_cv <- vfold_cv(bt_train,
                   v = 5,
                   strata = medv)

ข้อที่ 3. Grid ซึ่งเราจะใช้การจับคู่แบบสุ่มผ่าน grid_random():

# Set seed for reproducibility
set.seed(2025)

# Define the grid for tuning
hpt_grid <- grid_random(cost_complexity(range = c(-5, 0), trans = log10_trans()),
                        tree_depth(range = c(1, 20)),
                        min_n(range = c(2, 50)),
                        size = 20)

ข้อที่ 4. Metrics ซึ่งเราต้องเก็บไว้ใน metric_set():

# Define metrics
hpt_metrics = metric_set(mae,
                         rmse)

2️⃣ Tune

หลังจากเตรียม input ครบแล้ว เราสามารถ tune model ได้ด้วย tune_grid():

# Tune the model
dt_tune_results <- tune_grid(bt_wfl_tune,
                             resamples = hpt_cv,
                             grid = hpt_grid,
                             metrics = hpt_metrics)

3️⃣ Select

หลังจาก tune แล้ว เราสามารถดูค่า hyperparametres ที่ดีที่สุดตามตัวบ่งชี้ที่เราเลือกได้ด้วย show_best()

โดยในตัวอย่าง เราจะลองเลือก hyperparametres ที่ดีที่สุด 5 ชุดแรก โดยดูจาก RMSE:

# Show the best model
show_best(dt_tune_results,
          metric = "rmse",
          n = 5)

ผลลัพธ์:

# A tibble: 5 × 9
  cost_complexity tree_depth min_n .metric .estimator  mean     n
            <dbl>      <int> <int> <chr>   <chr>      <dbl> <int>
1       0.00311            3     5 rmse    standard    4.55     5
2       0.0167            16    12 rmse    standard    4.70     5
3       0.00150            5    28 rmse    standard    4.75     5
4       0.0000436         18    14 rmse    standard    4.88     5
5       0.0000345         11    39 rmse    standard    4.88     5
# ℹ 2 more variables: std_err <dbl>, .config <chr>

จากนั้น เราสามารถเลือกค่า hyperparametres ที่ดีที่สุดได้ด้วย select_best():

# Select the best model
dt_best_params <- select_best(dt_tune_results,
                              metric = "rmse")

4️⃣ Fit

ในขั้นที่ 4 เราจะใส่ค่า hyperparametres ที่เลือกมาเข้าไปใน workflow object ผ่าน finalize_workflow():

# Finalise the best workflow
dt_wkl_best <- finalize_workflow(bt_wfl_tune,
                                 dt_best_params)

จากนั้น train model ด้วย last_fit():

# Fit the best model
dt_best_fit <- last_fit(dt_wkl_best,
                        split = bt_split,
                        metrics = metric_set(mae, rmse))

5️⃣ Evaluate

สุดท้าย เราจะทดสอบ model ด้วย collect_predictions() และ collect_metrics():

# Collect predictions
predictions_best <- collect_predictions(dt_best_fit)

# Print predictions
predictions_best

ผลลัพธ์:

# A tibble: 103 × 5
   .pred id                .row  medv .config             
   <dbl> <chr>            <int> <dbl> <chr>               
 1  32.3 train/test split     3  34.7 Preprocessor1_Model1
 2  32.3 train/test split     4  33.4 Preprocessor1_Model1
 3  11.8 train/test split     9  16.5 Preprocessor1_Model1
 4  17.2 train/test split    25  15.6 Preprocessor1_Model1
 5  22.6 train/test split    40  30.8 Preprocessor1_Model1
 6  22.6 train/test split    53  25   Preprocessor1_Model1
 7  32.3 train/test split    56  35.4 Preprocessor1_Model1
 8  22.6 train/test split    79  21.2 Preprocessor1_Model1
 9  22.6 train/test split    82  23.9 Preprocessor1_Model1
10  34.5 train/test split    99  43.8 Preprocessor1_Model1
# ℹ 93 more rows
# ℹ Use `print(n = ...)` to see more rows

และ

# Collect metrics
metrics_best <- collect_metrics(dt_best_fit)

# Print metrics
metrics_best

ผลลัพธ์:

# A tibble: 2 × 4
  .metric .estimator .estimate .config             
  <chr>   <chr>          <dbl> <chr>               
1 mae     standard        3.48 Preprocessor1_Model1
2 rmse    standard        4.68 Preprocessor1_Model1

จะเห็นว่า hyperparametre tuning ทำให้ model ของเรามีประสิทธิภาพมากขึ้น เพราะมี MAE (3.70 vs 3.48) และ RMSE (5.13 vs 4.68) ที่ลดลง

😎 Summary

ในบทความนี้ เราได้ดูวิธีสร้าง ประเมิน และปรับ ML model ด้วย tidymodels ซึ่ง functions ต่าง ๆ ที่เราได้เรียนรู้สรุปตาม ML phase ได้ดังนี้:

Pre-processing:

Function	Description
`initial_split()`	แบ่งข้อมูล
`training()`	สร้าง training set
`testing()`	สร้าง test set
`recipe()`	กำหนดตัวแปรต้น ตัวแปรตาม
`step_*()`	กำหนดขั้นการแปลงข้อมูล
`prep()`	เตรียม recipe
`bake()`	เตรียมข้อมูล

Modelling:

Function	Description
`decision_tree()`	สร้าง decision tree
`set_engine()`	เรียกใช้ ML engine
`set_mode()`	กำหนดประเภท model
`fit()`	train model
`last_fit()`	train, ทำนาย, และประเมิน model

Post-processing:

Function	Description
`mae()`	คำนวณ MAE
`rmse()`	คำนวณ RMSE
`metric_set()`	กำหนดชุดตัวชี้วัด
`collect_predictions()`	เรียกดูคำทำนาย
`collect_metrics()`	เรียกดูตัวชี้วัด
`tune()`	กำหนด hyperparametres ที่ต้องการ tune
`vfold_cv()`	สร้าง k-fold cross-validation
`grid_random()`	สุ่มสร้างค่า hyperparametres ที่ต้องการทดสอบ
`tune_grid()`	tune model

All:

Function	Description
`workflow()`	รวม recipe และ model ไว้ด้วยกัน

📚 Further Reading

สำหรับคนที่สนใจ สามารถศึกษาเกี่ยวกับ tidymodels เพิ่มเติมได้ที่:

💪 Example Project

ดูตัวอย่างการใช้งาน tidymodels เพื่อสำหรับ data analytics project ได้ที่ Exploring & Predicting Employee Attrition With Machine Learning in R

😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub:

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-06-19

caret: แนะนำวิธีใช้แพ็กเกจ machine learning ที่ทรงพลังในภาษา R — ตัวอย่างการทำนายราคาบ้านใน BostonHousing dataset

caret เป็น package ยอดนิยมในภาษา R ในทำ machine learning (ML)

caret ย่อมาจาก Classification And REgression Training และเป็น package ที่ถูกออกแบบมาช่วยให้การทำ ML เป็นเรื่องง่ายโดยมี functions สำหรับทำงานกับ ML workflow เบ็ดเสร็จใน package เดียว

ในบทความนี้ เราจะไปดูตัวอย่างการใช้ caret กับ BostonHousing จาก mlbench กัน

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

🏠 Dataset: BostonHousing

BostonHousing เป็นชุดข้อมูลบ้านใน Boston, Massachusetts, US จาก mlbench package และประกอบด้วยข้อมูลบ้าน 14 columns ดังนี้:

No.	Column	Description
1	`crim`	ระดับอาชญากรรมในแต่ละเขต
2	`zn`	สัดส่วนพื้นที่อาศัย
3	`indus`	สัดส่วนธุรกิจที่เป็น non-retail ในแต่ละเขต
4	`chas`	เป็นพื้นที่ติดกับ Charles River ไหม (1 = ติด, 0 = ไม่ติด)
5	`nox`	ระดับ nitrogen oxide
6	`rm`	จำนวนห้องโดยเฉลี่ย
7	`age`	สัดส่วย unit ที่มีคนเข้าอยู่ ซึ่งถูกสร้างก่อนปี ค.ศ. 1940
8	`dis`	ระยะทางจากพื้นที่ทำงานในเมืองบอสตัน
9	`rad`	ระดับการเข้าถึง radial highways
10	`tax`	ภาษีโรงเรือน
11	`ptratio`	สัดส่วนนักเรียนต่อครูในแต่ละเขต
12	`black`	สัดส่วนผู้อยู่อาศัยที่เป็นคนผิวดำ
13	`lstat`	สัดส่วนของประชากรที่มีฐานะยากจน
14	`medv`	ราคากลางของบ้านที่มีผู้อยู่อาศัย

เป้าหมายในการทำงานกับ BostonHousing คือ สร้าง model เพื่อทำนายราคาบ้าน (medv)

⬇️ Load BostonHousing

เราสามารถโหลด BostonHousing เพื่อนำมาใช้งานได้แบบนี้:

# Install mlbench
install.packages("mlbench")

# Load mlbench
library(caret)
library(mlbench)

# Load BostonHousing
data("BostonHousing")

เราสามารถดูตัวอย่างข้อมูล BostonHousing ด้วย head():

# Preview
head(BostonHousing)

ผลลัพธ์:

     crim zn indus chas   nox    rm  age    dis rad tax ptratio      b lstat medv
1 0.00632 18  2.31    0 0.538 6.575 65.2 4.0900   1 296    15.3 396.90  4.98 24.0
2 0.02731  0  7.07    0 0.469 6.421 78.9 4.9671   2 242    17.8 396.90  9.14 21.6
3 0.02729  0  7.07    0 0.469 7.185 61.1 4.9671   2 242    17.8 392.83  4.03 34.7
4 0.03237  0  2.18    0 0.458 6.998 45.8 6.0622   3 222    18.7 394.63  2.94 33.4
5 0.06905  0  2.18    0 0.458 7.147 54.2 6.0622   3 222    18.7 396.90  5.33 36.2
6 0.02985  0  2.18    0 0.458 6.430 58.7 6.0622   3 222    18.7 394.12  5.21 28.7

และดูโครงสร้างชุดข้อมูลได้ด้วย str():

# View the structure
str(BostonHousing)

ผลลัพธ์:

'data.frame':	506 obs. of  14 variables:
 $ crim   : num  0.00632 0.02731 0.02729 0.03237 0.06905 ...
 $ zn     : num  18 0 0 0 0 0 12.5 12.5 12.5 12.5 ...
 $ indus  : num  2.31 7.07 7.07 2.18 2.18 2.18 7.87 7.87 7.87 7.87 ...
 $ chas   : Factor w/ 2 levels "0","1": 1 1 1 1 1 1 1 1 1 1 ...
 $ nox    : num  0.538 0.469 0.469 0.458 0.458 0.458 0.524 0.524 0.524 0.524 ...
 $ rm     : num  6.58 6.42 7.18 7 7.15 ...
 $ age    : num  65.2 78.9 61.1 45.8 54.2 58.7 66.6 96.1 100 85.9 ...
 $ dis    : num  4.09 4.97 4.97 6.06 6.06 ...
 $ rad    : num  1 2 2 3 3 3 5 5 5 5 ...
 $ tax    : num  296 242 242 222 222 222 311 311 311 311 ...
 $ ptratio: num  15.3 17.8 17.8 18.7 18.7 18.7 15.2 15.2 15.2 15.2 ...
 $ b      : num  397 397 393 395 397 ...
 $ lstat  : num  4.98 9.14 4.03 2.94 5.33 ...
 $ medv   : num  24 21.6 34.7 33.4 36.2 28.7 22.9 27.1 16.5 18.9 ...

🔧 Build a Predictive Model

การสร้าง model เพื่อทำนายราคาบ้านมีอยู่ 4 ขั้นตอน ได้แก่:

Split the data
Preprocess the data
Train the model
Evaluate the model

เราไปดูการใช้งาน caret สำหรับการทำงานแต่ละขั้นกัน

🪓 Step 1. Split the Data

ในขั้นแรก เราจะแบ่ง BostonHousing ออกเป็น 2 ชุด:

Training set สำหรับสร้าง model
Test set สำหรับประเมินประสิทธิภาพของ model

caret มี function สำหรับแบ่งชุดข้อมูล ซึ่งได้แก่ createDataPartition() ซึ่งต้องการ 3 arguments หลัก ดังนี้:

createDataPartition(y, p, list)

y = ตัวแปรตามที่เราต้องการทำนาย
p = สัดส่วนข้อมูลที่เราต้องการแบ่งให้กับ training set
list = ต้องการผลลัพธ์เป็น list (TRUE) หรือ matrix (FALSE)

สำหรับ BostonHousing เราจะแบ่ง 70% เป็น training set และ 30% เป็น test set แบบนี้:

# Set seed for reproducibility
set.seed(888)

# Get train index
train_index <- createDataPartition(BostonHousing$medv, # Specify the outcome
                                   p = 0.7, # Set aside 70% for training set
                                   list = FALSE) # Return as matrix

# Create training set
bt_train <- BostonHousing[train_index, ]

# Create test set
bt_test <- BostonHousing[-train_index, ]

เราสามารถดูจำนวนข้อมลใน training และ test sets ได้ด้วย nrow():

# Check the results
cat("Rows in training set:", nrow(bt_train), "\\n")
cat("Rows in test set:", nrow(bt_test))

ผลลัพธ์:

Rows in training set: 356
Rows in test set: 150

🍳 Step 2. Preprocess the Data

ในขั้นที่ 2 เราจะเตรียมข้อมูลเพื่อใช้ในการสร้าง model

caret มี function ที่ช่วยให้เราเตรียมข้อมูลได้ง่าย นั่นคือ preProcess() ซึ่งต้องการ 2 arguments หลัก ดังนี้:

preProcess(x, method)

x = ชุดข้อมูลที่เราต้องการใช้
method = วิธีการเตรียมข้อมูลที่เราต้องการ

สำหรับตัวอย่างนี้ เราจะใช้ 2 วิธีในการเตรียมข้อมูล ได้แก่:

centre หรือการปรับ mean ให้เท่ากับ 0
scale หรือการปรับ standard deviation ให้เท่ากับ 1

เราสามารถใช้ preProcess() ได้ดังนี้:

# Learn preprocessing on training set
ppc <- preProcess(bt_train[, -14], # Select all predictors
                  method = c("center", "scale"))  # Centre and scale

ตอนนี้ เราจะได้วิธีการเตรียมข้อมูลมาแล้ว ซึ่งเราจะต้องนำไปปรับใช้กับ training และ test sets ด้วย predict() และ cbind() แบบนี้:

Training set:

# Apply preprocessing to training set
bt_train_processed <- predict(ppc,
                              bt_train[, -14])

# Combine the preprocessed training set with outcome
bt_train_processed <- cbind(bt_train_processed,
                            medv = bt_train$medv)

Test set:

# Apply preprocessing to test set
bt_test_processed <- predict(ppc,
                             bt_test[, -14])

# Combine the preprocessed test set with outcome
bt_test_processed <- cbind(bt_test_processed,
                            medv = bt_test$medv)

ตอนนี้ training และ test sets ก็พร้อมที่จะใช้ในการสร้าง model แล้ว

👟 Step 3. Train the Model

ในขั้นที่ 3 เราจะสร้าง model กัน

โดยในตัวอย่าง เราจะลองสร้าง k-nearest neighbor (KNN) model ซึ่งทำนายข้อมูลด้วยการดูข้อมูลที่อยู่ใกล้เคียง

ในการสร้าง model, caret มี function ที่ทรงพลังและใช้งานง่ายอยู่ นั่นคือ train() ซึ่งต้องการ 5 arguments หลัก ดังนี้:

train(form, data, method, trControl, tuneGrid)

form = สูตรที่ระบุตัวแปรต้นและตัวแปรตาม (ตัวแปรตาม ~ ตัวแปรต้น)
data = ชุดข้อมูลที่จะใช้สร้าง model
method = engine ที่จะใช้ในการสร้าง model (ดูรายการ engine ทั้งหมดได้ที่ The caret Package — Available Models)
trControl = ค่าที่กำหนดการสร้าง model (ต้องใช้ function ชื่อ trainControl() ในการกำหนด)
tuneGrid = data frame ที่กำหนดค่า hyperparametre เพื่อทำ model tuning และหา model ที่ดีที่สุด

เราสามารถใช้ train() เพื่อสร้าง KNN model ในการทำนายราคาบ้านได้แบบนี้:

# Define training control:
# use k-fold cross-validation where k = 5
trc <- trainControl(method = "cv",
                    number = 5)

# Define grid:
# set k as odd numbers between 3 and 13
grid <- data.frame(k = seq(from = 3,
                           to = 13,
                           by = 2))

# Train the model
knn_model <- train(medv ~ ., # Specify the formula
                   data = bt_train_processed, # Use training set
                   method = "kknn", # Use knn engine
                   trControl = trc, # Specify training control
                   tuneGrid = grid) # Use grid to tune the model

เราสามารถดูรายละเอียดของ model ได้ดังนี้:

# Print the model
knn_model

ผลลัพธ์:

k-Nearest Neighbors 

356 samples
 13 predictor

No pre-processing
Resampling: Cross-Validated (5 fold) 
Summary of sample sizes: 284, 286, 286, 284, 284 
Resampling results across tuning parameters:

  k   RMSE      Rsquared   MAE     
   3  4.357333  0.7770080  2.840630
   5  4.438162  0.7760085  2.849984
   7  4.607954  0.7610468  2.941034
   9  4.683062  0.7577702  2.972661
  11  4.771317  0.7508908  3.043617
  13  4.815444  0.7524266  3.053415

RMSE was used to select the optimal model using the smallest value.
The final value used for the model was k = 3.

📈 Step 4. Evaluate the Model

ในขั้นสุดท้าย เราจะประเมินความสามารถของ model ในการทำนายราคาบ้านกัน

สำหรับ regression model, caret มี 3 functions ที่ช่วยในการประเมิน ได้แก่:

MAE() เพื่อหา mean absolute error (MAE) หรือค่าความคลาดเคลื่อนแบบสัมบูรณ์
RMSE() เพื่อหา root mean square error (RMSE) หรือค่าความคลาดเคลื่อนแบบกำลังสอง
R2() เพื่อหา R-squared หรือสัดส่วนของตัวแปรตามที่ model สามารถทำนายได้

เราสามารถใช้ทั้ง 3 functions ได้ดังนี้:

# Make predictions
predictions <- predict(knn_model,
                       newdata = bt_test_processed)

# Calculate MAE
mae <- MAE(predictions,
           bt_test_processed$medv)

# Calculate RMSE
rmse <- RMSE(predictions,
             bt_test_processed$medv)

# Calculate R squared
r2 <- R2(predictions,
         bt_test_processed$medv)

จากนั้น เราสามารถรวมผลลัพธ์ไว้ใน data frame เพื่อแสดงผลได้:

# Combine the results
results <- data.frame(Model = "KNN",
                      MAE = round(mae, 2),
                      RMSE = round(rmse, 2),
                      R_Squared = round(r2, 2))

# Print the results
results

ผลลัพธ์:

  Model  MAE RMSE R_Squared
1   KNN 2.67 3.74      0.86

จะเห็นได้ว่า model ของเรามีความคลาดเคลื่อนโดยเฉลี่ย 2.67 (MAE) และสามารถทำนาย 86% ของราคาบ้านทั้งหมดได้ (R_Squared)

💪 Summary

ในบทความนี้ เราได้ไปทำความรู้จักกับ caret ซึ่งเป็น package ในการสร้าง ML model ที่ทรงพลังและใช้งานง่าย โดยเราได้เรียนรู้ functions ในการทำงานดังนี้:

ML Workflow	Function	For
Split the data	`createDataPartition()`	สร้าง index สำหรับแบ่งข้อมูล
Preprocess the data	`preProcess()`	เตรียมข้อมูลสำหรับสร้าง model
Train the model	`train()`	สร้าง model
Evaluate the model	`MAE()`	หาค่า MAE
Evaluate the model	`RMSE()`	หาค่า RMSE
Evaluate the model	`R2()`	หาค่า R-squared

📚 Further Readings

สำหรับคนที่สนใจ สามารถอ่านเพิ่มเติมเกี่ยวกับ caret ได้ที่:

🐱 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-06-12

สอนปลูกต้นไม้ในภาษา R (ภาค 2): วิธีสร้าง ประเมิน และปรับทูน random forest model ด้วย ranger package – ตัวอย่างการทำนายระดับการกินน้ำมันของรถใน mpg dataset

ในบทความนี้ เราจะไปทำความรู้จักกับ random forest รวมทั้งการสร้าง ประเมิน และปรับทูน random forest model ด้วย ranger package ในภาษา R

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

🌲 Random Forest Model คืออะไร?

Random forest model เป็น tree-based model ซึ่งสุ่มสร้าง decision trees ขึ้นมาหลาย ๆ ต้น (forest) และใช้ผลลัพธ์ในภาพรวมเพื่อทำนายข้อมูลสุดท้าย:

Regression task: หาค่าเฉลี่ยของผลลัพธ์จากทุกต้น
Classification task: ดูผลลัพธ์ที่เป็นเสียงโหวตข้างมาก

Random forest เป็น model ที่ทรงพลัง เพราะใช้ผลรวมของหลาย ๆ decision trees แม้ว่า decision tree แต่ละต้นจะมีความสามารถในการทำนายนอยก็ตาม

💻 Random Forest Models ในภาษา R

ในภาษา R เรามี 2 packages ที่นิยมใช้สร้าง random forest model ได้แก่:

randomForest ซึ่งเป็น package ที่มีลูกเล่น แต่เก่ากว่า
ranger ซึ่งใหม่กว่า ประมวลผลได้เร็วกว่า และใช้งานง่ายกว่า

ในบทความก่อน เราดูวิธีการใช้ randomForest แล้ว

ในบทความนี้ เราจะไปดูวิธีใช้ ranger โดยใช้ mpg dataset เป็นตัวอย่างกัน

🚗 mpg Dataset

mpg dataset เป็น dataset จาก ggplots2 package และมีข้อมูลของรถ 38 รุ่น จากช่วงปี ช่วง ค.ศ. 1999 ถึง 2008 ทั้งหมด 11 columns ดังนี้:

No.	Column	Description
1	`manufacturer`	ผู้ผลิต
2	`model`	รุ่นรถ
3	`displ`	ขนาดถังน้ำมัน (ลิตร)
4	`year`	ปีที่ผลิต
5	`cyl`	จำนวนลูกสูบ
6	`trans`	ประเภทเกียร์
7	`drv`	ประเภทล้อขับเคลื่อน
8	`cty`	ระดับการกินน้ำมันเวลาวิ่งในเมือง
9	`hwy`	ระดับการกินน้ำมันเวลาวิ่งบน highway
10	`fl`	ประเภทน้ำมัน
11	`class`	ประเภทรถ

ในบทความนี้ เราจะลองใช้ ranger เพื่อทำนาย hwy กัน

เราสามารถเตรียม mpg เพื่อสร้าง random forest model ได้ดังนี้

โหลด dataset:

# Install ggplot2
install.packages("ggplot2")

# Load ggplot2
library(ggplot2)

# Load the dataset
data(mpg)

ดูตัวอย่าง dataset:

# Preview
head(mpg)

ผลลัพธ์:

# A tibble: 6 × 11
  manufacturer model displ  year   cyl trans      drv     cty   hwy fl    class  
  <chr>        <chr> <dbl> <int> <int> <chr>      <chr> <int> <int> <chr> <chr>  
1 audi         a4      1.8  1999     4 auto(l5)   f        18    29 p     compact
2 audi         a4      1.8  1999     4 manual(m5) f        21    29 p     compact
3 audi         a4      2    2008     4 manual(m6) f        20    31 p     compact
4 audi         a4      2    2008     4 auto(av)   f        21    30 p     compact
5 audi         a4      2.8  1999     6 auto(l5)   f        16    26 p     compact
6 audi         a4      2.8  1999     6 manual(m5) f        18    26 p     compact

สำรวจโครงสร้าง:

# View the structure
str(mpg)

ผลลัพธ์:

tibble [234 × 11] (S3: tbl_df/tbl/data.frame)
 $ manufacturer: chr [1:234] "audi" "audi" "audi" "audi" ...
 $ model       : chr [1:234] "a4" "a4" "a4" "a4" ...
 $ displ       : num [1:234] 1.8 1.8 2 2 2.8 2.8 3.1 1.8 1.8 2 ...
 $ year        : int [1:234] 1999 1999 2008 2008 1999 1999 2008 1999 1999 2008 ...
 $ cyl         : int [1:234] 4 4 4 4 6 6 6 4 4 4 ...
 $ trans       : chr [1:234] "auto(l5)" "manual(m5)" "manual(m6)" "auto(av)" ...
 $ drv         : chr [1:234] "f" "f" "f" "f" ...
 $ cty         : int [1:234] 18 21 20 21 16 18 18 18 16 20 ...
 $ hwy         : int [1:234] 29 29 31 30 26 26 27 26 25 28 ...
 $ fl          : chr [1:234] "p" "p" "p" "p" ...
 $ class       : chr [1:234] "compact" "compact" "compact" "compact" ...

จากผลลัพธ์จะเห็นว่า บาง columns (เช่น manufacturer, model) มีข้อมูลประเภท character ซึ่งเราควระเปลี่ยนเป็น factor เพื่อช่วยให้การสร้าง model มีประสิทธิภาพมากขึ้น:

# Convert character columns to factor

## Get character columns
chr_cols <- c("manufacturer", "model",
              "trans", "drv",
              "fl", "class")

## For-loop through the character columns
for (col in chr_cols) {
  mpg[[col]] <- as.factor(mpg[[col]])
}

## Check the results
str(mpg)

ผลลัพธ์:

tibble [234 × 11] (S3: tbl_df/tbl/data.frame)
 $ manufacturer: Factor w/ 15 levels "audi","chevrolet",..: 1 1 1 1 1 1 1 1 1 1 ...
 $ model       : Factor w/ 38 levels "4runner 4wd",..: 2 2 2 2 2 2 2 3 3 3 ...
 $ displ       : num [1:234] 1.8 1.8 2 2 2.8 2.8 3.1 1.8 1.8 2 ...
 $ year        : int [1:234] 1999 1999 2008 2008 1999 1999 2008 1999 1999 2008 ...
 $ cyl         : int [1:234] 4 4 4 4 6 6 6 4 4 4 ...
 $ trans       : Factor w/ 10 levels "auto(av)","auto(l3)",..: 4 9 10 1 4 9 1 9 4 10 ...
 $ drv         : Factor w/ 3 levels "4","f","r": 2 2 2 2 2 2 2 1 1 1 ...
 $ cty         : int [1:234] 18 21 20 21 16 18 18 18 16 20 ...
 $ hwy         : int [1:234] 29 29 31 30 26 26 27 26 25 28 ...
 $ fl          : Factor w/ 5 levels "c","d","e","p",..: 4 4 4 4 4 4 4 4 4 4 ...
 $ class       : Factor w/ 7 levels "2seater","compact",..: 2 2 2 2 2 2 2 2 2 2 ...

ตอนนี้ เราพร้อมที่จะนำ dataset ไปใช้งานกับ ranger แล้ว

🐣 ranger Basics

การใช้งาน ranger มีอยู่ 4 ขั้นตอน:

ติดตั้งและโหลด ranger
สร้าง training และ test sets
สร้าง random forest model
ทดสอบความสามารถของ model

1️⃣ ติดตั้งและโหลด ranger

ในครั้งแรกสุด ให้เราติดตั้ง ranger ด้วยคำสั่ง install.packages():

# Install
install.packages("ranger")

และทุกครั้งที่เราต้องการใช้งาน ranger ให้เราเรียกใช้งานด้วย library():

# Load
library(ranger)

2️⃣ สร้าง Training และ Test Sets

ในขั้นที่ 2 เราจะแบ่ง dataset เป็น 2 ส่วน ได้แก่:

Training set สำหรับสร้าง model (70% ของ dataset)
Test set สำหรับทดสอบ model (30% ของ dataset)

# Split the data

## Set seed for reproducibility
set.seed(123)

## Get training rows
train_rows <- sample(nrow(mpg),
                     nrow(mpg) * 0.7)

## Create a training set
train <- mpg[train_rows, ]

## Create a test set
test <- mpg[-train_rows, ]

3️⃣ สร้าง Random Forest Model

ในขั้นที่ 3 เราจะสร้าง random forest ด้วย ranger() ซึ่งต้องการ input หลัก 2 อย่าง ดังนี้:

ranger(formula, data)

formula: ระบุตัวแปรต้นและตัวแปรตามที่ใช้ในการสร้าง model
data: dataset ที่ใช้สร้าง model (เราจะใช้ training set กัน)

เราจะเรียกใช้ ranger() ดังนี้:

# Initial random forest model

## Set seed for reproducibility
set.seed(123)

## Train the model
rf_model <- ranger(hwy ~ .,
                   data = train)

Note: เราใช้ set.seed() เพื่อให้เราสามารถสร้าง model ซ้ำได้ เพราะ random forest มีการสร้าง decision trees แบบสุ่ม

เมื่อได้ model มาแล้ว เราสามารถดูรายละเอียดของ model ได้แบบนี้:

# Print the model
rf_model

ผลลัพธ์:

Ranger result

Call:
 ranger(hwy ~ ., data = train) 

Type:                             Regression 
Number of trees:                  500 
Sample size:                      163 
Number of independent variables:  10 
Mtry:                             3 
Target node size:                 5 
Variable importance mode:         none 
Splitrule:                        variance 
OOB prediction error (MSE):       1.682456 
R squared (OOB):                  0.9584596

ในผลลัพธ์ เราจะเห็นลักษณะต่าง ๆ ของ model เช่น ประเภท model (type) และ จำนวน decision trees ที่ถูกสร้างขึ้นมา (sample size)

4️⃣ ทดสอบความสามารถของ Model

สุดท้าย เราจะทดสอบความสามารถของ model ในการทำนายข้อมูล โดยเริ่มจากใช้ model ทำนายข้อมูลใน test set:

# Make predictions
preds <- predict(rf_model,
                 data = test)$predictions

จากนั้น คำนวณตัวบ่งชี้ความสามารถ (metric) ซึ่งสำหรับ regression model มีอยู่ 3 ตัว ได้แก่:

MAE (mean absolute error): ค่าเฉลี่ยความคลาดเคลื่อนแบบสัมบูรณ์ (ยิ่งน้อยยิ่งดี)
RMSE (root mean square error): ค่าเฉลี่ยความคาดเคลื่อนแบบยกกำลังสอง (ยิ่งน้อยยิ่งดี)
R squared: สัดส่วนข้อมูลที่อธิบายได้ด้วย model (ยิ่งมากยิ่งดี)

# Get errors
errors <- test$hwy - preds

# Calculate MAE
mae <- mean(abs(errors))

# Calculate RMSE
rmse <- sqrt(mean(errors^2))

# Calculate R squared
r_sq <- 1 - (sum((errors)^2) / sum((test$hwy - mean(test$hwy))^2))

# Print the results
cat("Initial model MAE:", round(mae, 2), "\n")
cat("Initial model RMSE:", round(rmse, 2), "\n")
cat("Initial model R squared:", round(r_sq, 2), "\n")

ผลลัพธ์:

Initial model MAE: 0.79
Initial model RMSE: 1.07
Initial model R squared: 0.95

⏲️ Hyperparametre Tuning

ranger มี hyperparametre มากมายที่เราสามารถปรับแต่งเพื่อเพิ่มประสิทธิภาพของ random forest model ได้ เช่น:

num.trees: จำนวน decision trees ที่จะสร้าง
mtry: จำนวนตัวแปรต้นที่จะถูกสุ่มไปใช้ในแต่ละ node
min.node.size: จำนวนข้อมูลขั้นต่ำที่แต่ละ node จะต้องมี

เราสามารถใช้ for loop เพื่อปรับหาค่า hyperparametre ที่ดีที่สุดได้ดังนี้:

# Define hyperparametres
ntree_vals <- c(300, 500, 700)
mtry_vals <- 2:5
min_node_vals <- c(1, 5, 10)

# Create a hyperparametre grid
grid <- expand.grid(num.trees = ntree_vals,
                    mtry = mtry_vals,
                    min.node.size = min_node_vals)

# Instantiate an empty data frame
hpt_results <- data.frame()

# For-loop through the hyperparametre grid
for (i in 1:nrow(grid)) {
  
  ## Get the combination
  params <- grid[i, ]
  
  ## Set seed for reproducibility
  set.seed(123)
  
  ## Fit the model
  model <- ranger(hwy ~ .,
                  data = train,
                  num.trees = params$num.trees,
                  mtry = params$mtry,
                  min.node.size = params$min.node.size)
  
  ## Make predictions
  preds <- predict(model,
                   data = test)$predictions
  
  ## Get errors
  errors <- test$hwy - preds
  
  ## Calculate MAE
  mae <- mean(abs(errors))
  
  ## Calculate RMSE
  rmse <- sqrt(mean(errors^2))
  
  ## Store the results
  hpt_results <- rbind(hpt_results,
                       cbind(params,
                             MAE = mae,
                             RMSE = rmse))
}

# View the results
hpt_results

ผลลัพธ์:

   num.trees mtry min.node.size       MAE      RMSE
1        300    2             1 0.8101026 1.0971836
2        500    2             1 0.8012484 1.0973957
3        700    2             1 0.8039271 1.1001252
4        300    3             1 0.7434543 1.0051344
5        500    3             1 0.7417985 1.0069989
6        700    3             1 0.7421666 1.0028184
7        300    4             1 0.6989314 0.9074216
8        500    4             1 0.7130704 0.9314843
9        700    4             1 0.7141147 0.9292718
10       300    5             1 0.7157657 0.9370918
11       500    5             1 0.7131899 0.9266787
12       700    5             1 0.7091556 0.9238312
13       300    2             5 0.8570125 1.1673637
14       500    2             5 0.8515116 1.1736009
15       700    2             5 0.8522571 1.1756648
16       300    3             5 0.7885005 1.0654548
17       500    3             5 0.7872713 1.0664734
18       700    3             5 0.7859149 1.0581331
19       300    4             5 0.7561500 0.9790160
20       500    4             5 0.7623437 0.9869463
21       700    4             5 0.7611660 0.9813048
22       300    5             5 0.7615190 0.9777769
23       500    5             5 0.7615861 0.9804616
24       700    5             5 0.7613151 0.9788333
25       300    2            10 0.9257704 1.2391377
26       500    2            10 0.9292344 1.2611164
27       700    2            10 0.9258555 1.2635794
28       300    3            10 0.8790601 1.1635695
29       500    3            10 0.8704461 1.1594165
30       700    3            10 0.8704562 1.1507016
31       300    4            10 0.8609516 1.0887466
32       500    4            10 0.8672105 1.0962367
33       700    4            10 0.8624934 1.0875710
34       300    5            10 0.8558867 1.0811168
35       500    5            10 0.8567463 1.0783473
36       700    5            10 0.8536824 1.0751511

จะเห็นว่า เราจะได้ MSE และ RMSE ของส่วนผสมระหว่างแต่ละ hyperparametre มา

เราสามารถใช้ ggplot() เพื่อช่วยเลือก hyperparametres ที่ดีที่สุดได้ดังนี้:

# Visualise the results
ggplot(hpt_results,
       aes(x = mtry,
           y = RMSE,
           color = factor(num.trees))) +
  
  ## Use scatter plot
  geom_point(aes(size = min.node.size)) +
  
  ## Set theme to minimal
  theme_minimal() +
  
  ## Add title, labels, and legends
  labs(title = "Hyperparametre Tuning Results",
       x = "mtry",
       y = "RMSE",
       color = "num.trees",
       size = "min.node.size")

ผลลัพธ์:

จากกราฟ จะเห็นได้ว่า hyperparametres ที่ดีที่สุด (มี RMSE น้อยที่สุด) คือ:

num.trees = 300
mtry = 4
min.node.size = 2.5

เมื่อได้ค่า hyperparametres แล้ว เราสามารถใส่ค่าเหล่านี้กลับเข้าไปใน model และทดสอบความสามารถได้เลย

สร้าง model:

# Define the best hyperparametres
best_num.tree <- 300
best_mtry <- 4
best_min.node.size <- 2.5

# Fit the model
rf_model_new <- ranger(hwy ~ .,
                       data = train,
                       num.tree = best_num.tree,
                       mtry = best_mtry,
                       min.node.size = best_min.node.size)

ทดสอบความสามารถ:

# Evaluate the model

## Make predictions
preds_new <- predict(rf_model_new,
                     data = test)$predictions

## Get errors
errors_new <- test$hwy - preds_new

## Calculate MAE
mae_new <- mean(abs(errors_new))

## Calculate RMSE
rmse_new <- sqrt(mean(errors_new^2))

## Calculate R squared
r_sq_new <- 1 - (sum((errors_new)^2) / sum((test$hwy - mean(test$hwy))^2))

## Print the results
cat("Final model MAE:", round(mae_new, 2), "\n")
cat("Final model RMSE:", round(rmse_new, 2), "\n")
cat("Final model R squared:", round(r_sq_new, 2), "\n")

ผลลัพธ์:

Final model MAE: 0.71
Final model RMSE: 0.93
Final model R squared: 0.96

เราสามารถเปรียบเทียบความสามารถของ model ล่าสุด (final model) กับ model ก่อนหน้านี้ (initial model) ได้:

# Compare the two models
model_comp <- data.frame(Model = c("Initial", "Final"),
                         MAE = c(round(mae, 2), round(mae_new, 2)),
                         RMSE = c(round(rmse, 2), round(rmse_new, 2)),
                         R_Squared = c(round(r_sq, 2), round(r_sq_new, 2)))

# Print
model_comp

ผลลัพธ์:

    Model  MAE RMSE R_Squared
1 Initial 0.85 1.08      0.95
2   Final 0.71 0.93      0.96

ซึ่งจะเห็นว่า model ใหม่สามารถทำนายข้อมูลได้ดีขึ้น เพราะมี MAE และ RMSE ที่ลดลง รวมทั้ง R squared ที่เพิ่มขึ้น

🍩 Bonus: Variable Importance

ส่งท้าย ในกรณีที่เราต้องการดูว่า ตัวแปรต้นไหนมีความสำคัญต่อการทำนายมากที่สุด เราสามารถใช้ importance argument ใน ranger() คู่กับ vip() จาก vip package ได้แบบนี้:

# Fit the model with importance
rf_model_new <- ranger(hwy ~ .,
                       data = train,
                       num.tree = best_num.tree,
                       mtry = best_mtry,
                       min.node.size = best_min.node.size,
                       importance = "permutation") # Add importance

# Install vip package
install.packages("vip")

# Load vip package
library(vip)

# Get variabe importance
vip(rf_model_new)  +
  
  ## Add title and labels
  labs(title = "Variable Importance - Final Random Forest Model",
       x = "Variables",
       y = "Importance") +
  
  ## Set theme to minimal
  theme_minimal()

ผลลัพธ์:

จากกราฟ จะเห็นได้ว่า ตัวแปรต้นที่สำคัญที่สุด 3 ตัว ได้แก่:

cty: ระดับการกินน้ำมันเวลาวิ่งในเมือง
displ: ขนาดถังน้ำมัน (ลิตร)
cyl: จำนวนลูกสูบ

😎 Summary

ในบทความนี้ เราได้ดูวิธีการใช้ ranger package เพื่อ:

สร้าง random forest model
ปรับทูน model

พร้อมวิธีการประเมิน model ด้วย predict() และการคำนวณ MAE, RMSE, และ R squared รวมทั้งดูความสำคัญของตัวแปรต้นด้วย vip package

😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-06-05

วิธีสร้าง linear regression ด้วย lm() ในภาษา R — ตัวอย่างการทำนายราคาเพชรใน diamonds dataset

Linear regression เป็นวิธีการทำนายข้อมูลด้วยสมการเส้นตรง:

y = a + bx

y = ตัวแปรตาม หรือข้อมูลที่ต้องการทำนาย
a = จุดตัดระหว่าง x และ y (intercept)
b = ค่าความชัด (slope)
x = ตัวแปรต้น

เนื่องจากเป็นเทคนิคที่ใช้งานและทำความเข้าใจได้ง่าย linear regression จึงเป็นวิธีที่นิยมใช้ในการทำนายข้อมูลในบริบทต่าง ๆ เช่น:

ทำนาย	จาก
กำไร	ค่าโฆษณา
ความสามารถของนักกีฬา	ชั่วโมงฝึกซ้อม
ความดันเลือด	ปริมาณยา + อายุ
ผลลิตทางการเกษตร	ปริมาณน้ำ + ปุ๋ย

ในบทความนี้ เราจะมาดูวิธีใช้ linear regression ในภาษา R กัน

ถ้าพร้อมแล้ว ไปเริ่มกันเลย

💎 Example Dataset: diamonds

ในบทความนี้ เราจะใช้ diamonds dataset เป็นตัวอย่างในการใช้ linear regression กัน

diamonds dataset เป็น built-in dataset จาก ggplot2 package ซึ่งมีข้อมูลเพชรมากกว่า 50,000 ตัวอย่าง และประกอบด้วย 10 columns ดังนี้:

No.	Column	Description
1	`price`	ราคา (ดอลล่าร์สหรัฐฯ)
2	`caret`	น้ำหนัก
3	`cut`	คุณภาพ
4	`color`	สี
5	`clarity`	ความใสของเพชร
6	`x`	ความยาว
7	`y`	ความกว้าง
8	`z`	ความลึก
9	`depth`	สัดส่วนความลึก
10	`table`	สัดส่วนความกว้างของยอดเพชรต่อส่วนที่กว้างที่สุด

เป้าหมายของเรา คือ ทำนายราคาเพชร (price)

⬇️ Load diamonds

ในการใช้งาน diamonds เราสามารถเรียกใช้งาน dataset ได้ดังนี้:

ขั้นที่ 1. ติดตั้งและโหลด ggplot2:

# Install
install.packages("ggplot2")

# Load
library(ggplot2)

ขั้นที่ 2. โหลด diamonds dataset:

# Load dataset
data(diamonds)

ขั้นที่ 3. ดูตัวอย่างข้อมูล 10 rows แรกใน dataset:

# Preview the dataset
head(diamonds, 10)

ผลลัพธ์:

# A tibble: 10 × 10
   carat cut       color clarity depth table price     x     y     z
   <dbl> <ord>     <ord> <ord>   <dbl> <dbl> <int> <dbl> <dbl> <dbl>
 1  0.23 Ideal     E     SI2      61.5    55   326  3.95  3.98  2.43
 2  0.21 Premium   E     SI1      59.8    61   326  3.89  3.84  2.31
 3  0.23 Good      E     VS1      56.9    65   327  4.05  4.07  2.31
 4  0.29 Premium   I     VS2      62.4    58   334  4.2   4.23  2.63
 5  0.31 Good      J     SI2      63.3    58   335  4.34  4.35  2.75
 6  0.24 Very Good J     VVS2     62.8    57   336  3.94  3.96  2.48
 7  0.24 Very Good I     VVS1     62.3    57   336  3.95  3.98  2.47
 8  0.26 Very Good H     SI1      61.9    55   337  4.07  4.11  2.53
 9  0.22 Fair      E     VS2      65.1    61   337  3.87  3.78  2.49
10  0.23 Very Good H     VS1      59.4    61   338  4     4.05  2.39

🍳 Prepare the Dataset

ก่อนจะทำนายราคาเพชรด้วย linear regression เราจะเตรียม diamonds dataset ใน 3 ขั้นตอนก่อน ได้แก่:

One-hot encoding
Log transformation
Split data

🪆 Step 1. One-Hot Encoding

ในกรณีที่ตัวแปรต้นที่เป็น categorical เราจะต้องแปลงตัวแปรเหล่านี้ให้เป็น numeric ก่อน ซึ่งเราสามารถทำได้ด้วย one-hot encoding ดังตัวอย่าง:

ก่อน one-hot encoding:

Data	Cut
1	Ideal
2	Good
3	Fair

หลัง one-hot encoding:

Data	Cut_Ideal	Cut_Good	Cut_Fair
1	1	0	0
2	0	1	0
3	0	0	1

ในภาษา R เราสามารถทำ one-hot encoding ได้ด้วย model.matrix() ดังนี้:

# Set option for one-hot encoding
options(contrasts = c("contr.treatment",
                      "contr.treatment"))

# One-hot encode
cat_dum <- model.matrix(~ cut + color + clarity - 1,
                        data = diamonds)

จากนั้น เราจะนำผลลัพธ์ที่ได้ไปรวมกับตัวแปรตามและตัวแปรต้นที่เป็น numeric:

# Combine one-hot-encoded categorical and numeric variables
dm <- cbind(diamonds[, c("carat",
                         "depth",
                         "table",
                         "x",
                         "y",
                         "z")],
            cat_dum,
            price = diamonds$price)

เราสามารถเช็กผลลัพธ์ของ one-hot encoding ได้ด้วย str():

# Check the results
str(dm)

ผลลัพธ์:

'data.frame':	53940 obs. of  25 variables:
 $ carat       : num  0.23 0.21 0.23 0.29 0.31 0.24 0.24 0.26 0.22 0.23 ...
 $ depth       : num  61.5 59.8 56.9 62.4 63.3 62.8 62.3 61.9 65.1 59.4 ...
 $ table       : num  55 61 65 58 58 57 57 55 61 61 ...
 $ x           : num  3.95 3.89 4.05 4.2 4.34 3.94 3.95 4.07 3.87 4 ...
 $ y           : num  3.98 3.84 4.07 4.23 4.35 3.96 3.98 4.11 3.78 4.05 ...
 $ z           : num  2.43 2.31 2.31 2.63 2.75 2.48 2.47 2.53 2.49 2.39 ...
 $ cutFair     : num  0 0 0 0 0 0 0 0 1 0 ...
 $ cutGood     : num  0 0 1 0 1 0 0 0 0 0 ...
 $ cutVery Good: num  0 0 0 0 0 1 1 1 0 1 ...
 $ cutPremium  : num  0 1 0 1 0 0 0 0 0 0 ...
 $ cutIdeal    : num  1 0 0 0 0 0 0 0 0 0 ...
 $ colorE      : num  1 1 1 0 0 0 0 0 1 0 ...
 $ colorF      : num  0 0 0 0 0 0 0 0 0 0 ...
 $ colorG      : num  0 0 0 0 0 0 0 0 0 0 ...
 $ colorH      : num  0 0 0 0 0 0 0 1 0 1 ...
 $ colorI      : num  0 0 0 1 0 0 1 0 0 0 ...
 $ colorJ      : num  0 0 0 0 1 1 0 0 0 0 ...
 $ claritySI2  : num  1 0 0 0 1 0 0 0 0 0 ...
 $ claritySI1  : num  0 1 0 0 0 0 0 1 0 0 ...
 $ clarityVS2  : num  0 0 0 1 0 0 0 0 1 0 ...
 $ clarityVS1  : num  0 0 1 0 0 0 0 0 0 1 ...
 $ clarityVVS2 : num  0 0 0 0 0 1 0 0 0 0 ...
 $ clarityVVS1 : num  0 0 0 0 0 0 1 0 0 0 ...
 $ clarityIF   : num  0 0 0 0 0 0 0 0 0 0 ...
 $ price_log   : num  5.79 5.79 5.79 5.81 5.81 ...

ตอนนี้ ตัวแปรต้นที่เป็น categorical ถูกแปลงเป็น numeric ทั้งหมดแล้ว

📈 Step 2. Log Transformation

ในกรณีที่ตัวแปรตามมีการกระจายตัว (distribution) ไม่ปกติ linear regression ทำนายข้อมูลได้ไม่เต็มประสิทธิภาพนัก

เราสามารถตรวจสอบการกระจายตัวของตัวแปรตามได้ด้วย ggplot():

# Check the distribution of `price`
ggplot(dm,
       aes(x = price)) +
  
  ## Instantiate a histogram
  geom_histogram(binwidth = 100,
                 fill = "skyblue3") +
  
  ## Add text elements
  labs(title = "Distribution of Price",
       x = "Price",
       y = "Count") +
  
  ## Set theme to minimal
  theme_minimal()

ผลลัพธ์:

จากกราฟ เราจะเห็นได้ว่า ตัวแปรตามมีการกระจายตัวแบบเบ้ขวา (right-skewed)

ดังนั้น ก่อนจะใช้ linear regression เราจะต้องแปรตัวแปรตามให้มีการกระจายตัวแบบปกติ (normal distribution) ก่อน ซึ่งเราสามารถทำได้ด้วย log transformation ดังนี้:

# Log-transform `price`
dm$price_log <- log(dm$price)

# Drop `price`
dm$price <- NULL

หลัง log transformation เราสามารถเช็กการกระจายตัวด้วย ggplot() อีกครั้ง:

# Check the distribution of logged `price`
ggplot(dm,
       aes(x = price_log)) +
  
  ## Instantiate a histogram
  geom_histogram(fill = "skyblue3") +
  
  ## Add text elements
  labs(title = "Distribution of Price After Log Transformation",
       x = "Price (Logged)",
       y = "Count") +
  
  ## Set theme to minimal
  theme_minimal()

ผลลัพธ์:

จะเห็นได้ว่า การกระจายตัวของตัวแปรตามใกล้เคียงกับการกระจายตัวแบบปกติมากขึ้นแล้ว

🚄 Step 3. Split the Data

ในขั้นสุดท้ายก่อนใช้ linear regression เราจะแบ่งข้อมูลออกเป็น 2 ชุด:

Training set สำหรับสร้าง linear regression model
Test set สำหรับประเมินความสามารถของ linear regression model

ในบทความนี้ เราจะแบ่ง 80% ของ dataset เป็น training set และ 20% เป็น test set:

# Split the data

## Set seed for reproducibility
set.seed(181)

## Training index
train_index <- sample(nrow(dm),
                      0.8 * nrow(dm))

## Create training set
train_set <- dm[train_index, ]

## Create test set
test_set <- dm[-train_index, ]

ตอนนี้ เราพร้อมที่จะสร้าง linear regression model กันแล้ว

🏷️ Linear Regression Modelling

การสร้าง linear regression model มีอยู่ 3 ขั้นตอน ได้แก่:

Fit the model
Make predictions
Evaluate the model performance

💪 Step 1. Fit the Model

ในขั้นแรก เราจะสร้าง model ด้วย lm() ซึ่งต้องการ input 2 อย่าง:

lm(formula, data)

formula = สูตรการทำนาย โดยเราต้องกำหนดตัวแปรต้นและตัวแปรตาม
data = ชุดข้อมูลที่ใช้สร้าง model

ในการทำนายราคาเพชร เราจะใช้ lm() แบบนี้:

# Fit the model
linear_reg <- lm(price_log ~ .,
                 data = train_set)

อธิบาย code:

price_log ~ . หมายถึง ทำนายราคา (price_log) ด้วยตัวแปรต้นทั้งหมด (.)
data = train_set หมายถึง เรากำหนดชุดข้อมูลที่ใช้เป็น training set

เราสามารถดูข้อมูลของ model ได้ด้วย summary():

# View the model
summary(linear_reg)

ผลลัพธ์:

Call:
lm(formula = price_log ~ ., data = train_set)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.2093 -0.0930  0.0019  0.0916  9.8935 

Coefficients: (1 not defined because of singularities)
                 Estimate Std. Error  t value Pr(>|t|)    
(Intercept)    -2.7959573  0.0705854  -39.611  < 2e-16 ***
carat          -0.5270039  0.0086582  -60.867  < 2e-16 ***
depth           0.0512357  0.0008077   63.437  < 2e-16 ***
table           0.0090154  0.0005249   17.175  < 2e-16 ***
x               1.1374016  0.0055578  204.651  < 2e-16 ***
y               0.0290584  0.0031345    9.271  < 2e-16 ***
z               0.0340298  0.0054896    6.199 5.73e-10 ***
cutFair        -0.1528658  0.0060005  -25.476  < 2e-16 ***
cutGood        -0.0639105  0.0036547  -17.487  < 2e-16 ***
`cutVery Good` -0.0313800  0.0025724  -12.199  < 2e-16 ***
cutPremium     -0.0451760  0.0026362  -17.137  < 2e-16 ***
cutIdeal               NA         NA       NA       NA    
colorE         -0.0573940  0.0032281  -17.779  < 2e-16 ***
colorF         -0.0892633  0.0032654  -27.336  < 2e-16 ***
colorG         -0.1573861  0.0032031  -49.136  < 2e-16 ***
colorH         -0.2592763  0.0034037  -76.175  < 2e-16 ***
colorI         -0.3864526  0.0038360 -100.742  < 2e-16 ***
colorJ         -0.5258789  0.0047183 -111.455  < 2e-16 ***
claritySI2      0.4431577  0.0079170   55.976  < 2e-16 ***
claritySI1      0.6087513  0.0078819   77.234  < 2e-16 ***
clarityVS2      0.7523161  0.0079211   94.976  < 2e-16 ***
clarityVS1      0.8200656  0.0080463  101.918  < 2e-16 ***
clarityVVS2     0.9381319  0.0082836  113.252  < 2e-16 ***
clarityVVS1     1.0033931  0.0085098  117.910  < 2e-16 ***
clarityIF       1.0898015  0.0092139  118.277  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.1825 on 43128 degrees of freedom
Multiple R-squared:  0.9677,	Adjusted R-squared:  0.9676 
F-statistic: 5.611e+04 on 23 and 43128 DF,  p-value: < 2.2e-16

Note: ดูวิธีการอ่านผลลัพธ์ได้ที่ Explaining the lm() Summary in R และ Understanding Linear Regression Output in R

🔮 Step 2. Make Predictions

ในขั้นที่สอง เราจะใช้ model เพื่อทำนายราคาด้วย predict():

# Predict in the outcome space
pred <- exp(pred_log)

# Preview predictions
head(pred_log)

ผลลัพธ์:

       2        5        9       16       19       22 
5.828071 5.816460 6.111859 5.777434 5.865820 6.088356

จะเห็นว่า ราคาที่ทำนายยังอยู่ในรูป log ซึ่งเราต้องแปลงกลับเป็นราคาปกติด้วย exp():

# Predict in the outcome space
pred <- exp(pred_log)

# Preview predictions
head(pred)

ผลลัพธ์:

       2        5        9       16       19       22 
339.7028 335.7812 451.1766 322.9295 352.7713 440.6961

เราสามารถเปรียบเทียบราคาจริงกับราคาที่ทำนาย พร้อมความคลาดเคลื่อน ได้ดังนี้:

# Compare predictions to actual
results <- data.frame(actual = round(exp(test_set$price_log), 2),
                      predicted = round(pred, 2),
                      diff = round(exp(test_set$price_log) - pred, 2))

# Print results
head(results)

ผลลัพธ์:

   actual predicted    diff
2     326    339.70  -13.70
5     335    335.78   -0.78
9     337    451.18 -114.18
16    345    322.93   22.07
19    351    352.77   -1.77
22    352    440.70  -88.70

🎯 Step 3. Evaluate the Model Performance

ในขั้นสุดท้าย เราจะประเมิน model โดยใช้ 2 ตัวชี้วัด ได้แก่:

Mean absolute error (MAE): ค่าเฉลี่ยความคลาดเคลื่อนโดยสัมบูรณ์
Root mean squared error (RMSE): ค่าเฉลี่ยความคลาดเคลื่อนแบบยกกำลังสอง

ทั้งสองตัวคำนวณความแตกต่างระหว่างสิ่งที่ทำนายและข้อมูลจริง ยิ่ง MAE และ RMSE สูง ก็หมายความว่า การทำนายมีความคาดเคลื่อนมาก แสดงว่า model ทำงานได้ไม่ดีนัก

ในทางกลับกัน ถ้า MAE และ RMSE น้อย ก็แสดงว่า การทำนายใกล้เคียงกับข้อมูลจริง และ model มีความแม่นยำสูง

(Note: เรียนรู้ความแตกต่างระหว่าง MAE และ RMSE ได้ที่ Loss Functions in Machine Learning Explained)

เราสามารถคำนวณ MAE และ RMSE ได้ดังนี้:

# Calculate MAE
mae <- mean(abs(results$diff))

# Calculate RMSE
rmse <- sqrt(mean((results$diff)^2))

# Print the results
cat("MAE:", round(mae, 2), "\n")
cat("RMSE:", round(rmse, 2))

ผลลัพธ์:

MAE: 491.71
RMSE: 1123.68

จากผลลัพธ์ เราจะเห็นว่า โดยเฉลี่ย model ทำนายราคาคลาดเคลื่อนไปประมาณ 492 ดอลล่าร์ (MAE)

😎 Summary

ในบทความนี้ เราได้ดูวิธีการทำ linear regression ในภาษา R กัน

เราดูวิธีการเตรียมข้อมูลสำหรับ linear regression:

One-hot encoding ด้วย model.matrix()
Log transformation ด้วย log()
Split data ด้วย sample()

สร้าง linear regression model ด้วย lm() พร้อมประเมิน model ด้วย predict() และการคำนวณค่า MAE และ RMSE

😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

Correlation
t-tests
ANOVA
Reliability
Factor analysis

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists

2025-05-29

วิธีวิเคราะห์และแปลผล principal component analysis (PCA) ในภาษา R — ตัวอย่างการใช้ prcomp() เพื่อลดตัวแปรใน wine dataset
ในบทความนี้ เราจะมาทำความรู้จักกับ principal component analysis (PCA) ในภาษา R กัน
🧐 PCA คืออะไร?

PCA เป็น machine learning algorithm ประเภท unsupervised learning สำหรับลดตัวแปรในข้อมูล (dimensionality reduction) ในขณะที่ยังเก็บข้อมูลที่มีความสำคัญเอาไว้

ยกตัวอย่างเช่น เรามีข้อมูลลูกค้าที่มี 50 ตัวแปร เช่น อายุ เงินเดือน ประวัติการใช้จ่าย ประวัติการเป็นสมาชิก … เราสามารถใช้ PCA เพื่อลดจำนวนตัวแปรลง เช่น ลดเหลือ 5 ตัวแปร ที่ยังให้ข้อมูลเกี่ยวกับลูกค้าได้เทียบเท่ากับ 50 ตัวแปร

ในการทำงาน เราสามารถใช้ PCA เพื่อ:
- ลดเวลาและ resource ในการประมวลผล เพราะมีตัวแปรที่ต้องประมวลผลน้อยลง
- ช่วยทำความเข้าใจข้อมูล เพราะมีตัวแปรที่ต้องพิจารณาน้อยลง
- ใช้เตรียมข้อมูลสำหรับสร้าง machine learning models อื่น ๆ เช่น regression model, hierarchical regression
💻 PCA ในภาษา R

.

🧑‍💻 prcomp()

ในภาษา R เราสามารถทำ PCA ได้ด้วย prcomp() ซึ่งเป็น function ใน base R (ซึ่งหมายความว่า เราไม่ต้องติดตั้ง package เพิ่มเติม)

prcomp() ต้องการ 3 arguments ได้แก่:
```
prcomp(x, center, scale.)
```
- x คือ dataset ที่ต้องการใช้
- center คือ ตัวเลือกว่า เราจะลบ mean ออกจากข้อมูลดิบ เพื่อให้ dataset มีค่า mean เป็น 0 ไหม (recommend ให้เป็น TRUE)
- scale. คือ ตัวเลือกว่า เราจะหารข้อมูลดิบด้วย variance เพื่อให้ทุก column อยู่ในช่วงข้อมูลเดียวกันไหม (recommend ให้เป็น TRUE)
.

🍷 wine dataset

เราลองมาดูตัวอย่าง PCA กับ wine dataset จาก rattle package กัน

wine เป็นชุดข้อมูลที่มีลักษณะต่าง ๆ ของไวน์ เช่น ระดับแอลกอฮอล์ สี และความเข้มข้น

(Note: ดูข้อมูลเพิ่มเติมเกี่ยวกับ wine dataset ได้ที่ wine: The wine dataset from the UCI Machine Learning Repository.)

เราสามารถเรียกใช้ wine dataset ได้จาก rattle package ดังนี้:
```
# Install and load the package

## Install
install.packages("rattle")

## Load
library(rattle)

# -----------------------------------

# Load the dataset

## Load
data(wine)
```
เรียกดูตัวอย่างข้อมูลใน wine dataset:
```
# Preview
head(wine)
```
ผลลัพธ์:
```
  Type Alcohol Malic  Ash Alcalinity Magnesium Phenols Flavanoids Nonflavanoids Proanthocyanins Color  Hue Dilution Proline
1    1   14.23  1.71 2.43       15.6       127    2.80       3.06          0.28            2.29  5.64 1.04     3.92    1065
2    1   13.20  1.78 2.14       11.2       100    2.65       2.76          0.26            1.28  4.38 1.05     3.40    1050
3    1   13.16  2.36 2.67       18.6       101    2.80       3.24          0.30            2.81  5.68 1.03     3.17    1185
4    1   14.37  1.95 2.50       16.8       113    3.85       3.49          0.24            2.18  7.80 0.86     3.45    1480
5    1   13.24  2.59 2.87       21.0       118    2.80       2.69          0.39            1.82  4.32 1.04     2.93     735
6    1   14.20  1.76 2.45       15.2       112    3.27       3.39          0.34            1.97  6.75 1.05     2.85    1450
```
ดูโครงสร้างของ wine dataset:
```
# View the structure
str(wine)
```
ผลลัพธ์:
```
'data.frame':	178 obs. of  14 variables:
 $ Type           : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
 $ Alcohol        : num  14.2 13.2 13.2 14.4 13.2 ...
 $ Malic          : num  1.71 1.78 2.36 1.95 2.59 1.76 1.87 2.15 1.64 1.35 ...
 $ Ash            : num  2.43 2.14 2.67 2.5 2.87 2.45 2.45 2.61 2.17 2.27 ...
 $ Alcalinity     : num  15.6 11.2 18.6 16.8 21 15.2 14.6 17.6 14 16 ...
 $ Magnesium      : int  127 100 101 113 118 112 96 121 97 98 ...
 $ Phenols        : num  2.8 2.65 2.8 3.85 2.8 3.27 2.5 2.6 2.8 2.98 ...
 $ Flavanoids     : num  3.06 2.76 3.24 3.49 2.69 3.39 2.52 2.51 2.98 3.15 ...
 $ Nonflavanoids  : num  0.28 0.26 0.3 0.24 0.39 0.34 0.3 0.31 0.29 0.22 ...
 $ Proanthocyanins: num  2.29 1.28 2.81 2.18 1.82 1.97 1.98 1.25 1.98 1.85 ...
 $ Color          : num  5.64 4.38 5.68 7.8 4.32 6.75 5.25 5.05 5.2 7.22 ...
 $ Hue            : num  1.04 1.05 1.03 0.86 1.04 1.05 1.02 1.06 1.08 1.01 ...
 $ Dilution       : num  3.92 3.4 3.17 3.45 2.93 2.85 3.58 3.58 2.85 3.55 ...
 $ Proline        : int  1065 1050 1185 1480 735 1450 1290 1295 1045 1045 ...
```
Note: จะเห็นว่า ทุก column เป็นข้อมูลประเภท numeric ยกเว้น Type ที่ะเป็น factor

.

⚒️ PCA

หลังโหลด dataset แล้ว เราสามารถใช้ prcomp() เพื่อทำ PCA ได้ดังนี้:
```
# PCA
pca <- prcomp(wine[, -1],
              center = TRUE,
              scale. = TRUE)
```
อธิบาย code:
- เราใส่ wine[, -1] เพราะ PCA ใช้งานได้กับ column ที่เป็น numeric เท่านั้น เราเลย subset ข้อมูลเพื่อนำ Type ที่เป็น factor ออก
- เรากำหนด center = TRUE เพื่อให้ mean เท่ากับ 0
- กำหนดให้ scale. = TRUE เพื่อให้ข้อมูลอยู่ในช่วงข้อมูลเดียวกัน และป้องกันไม่ให้ข้อมูลที่อยู่ในมีช่วงข้อมูลกว้าง (เช่น ช่วง 1-100) มีผลต่อการวิเคราะห์มากกว่าข้อมูลที่มีช่วงแคบ (เช่น 1-10)
.

📈 ดูผลลัพธ์

เราสามารถดูผลลัพธ์ของ PCA ได้ 2 วิธี ได้แก่:
1. ดูค่าทางสถิติ
2. สร้างกราฟ
วิธีที่ 1. เราสามารถดูค่าทางสถิติได้ด้วย summary():
```
# Print the results
summary(pca)
```
ผลลัพธ์:
```
Importance of components:
                         PC1    PC2    PC3     PC4     PC5     PC6     PC7     PC8     PC9   PC10    PC11    PC12    PC13
Standard deviation     2.169 1.5802 1.2025 0.95863 0.92370 0.80103 0.74231 0.59034 0.53748 0.5009 0.47517 0.41082 0.32152
Proportion of Variance 0.362 0.1921 0.1112 0.07069 0.06563 0.04936 0.04239 0.02681 0.02222 0.0193 0.01737 0.01298 0.00795
Cumulative Proportion  0.362 0.5541 0.6653 0.73599 0.80162 0.85098 0.89337 0.92018 0.94240 0.9617 0.97907 0.99205 1.00000
```
ในผลลัพธ์ เราจะเห็นรายละเอียดดังนี้:
- จำนวนตัวแปรที่ PCA สร้างให้ (PC ซึ่งย่อมาจาก principal component) ซึ่งในตัวอย่างมีทั้งหมด 13 ตัวแปร
- Standard deviation (SD) ของแต่ละ PC
- สัดส่วนข้อมูลที่อธิบายได้ด้วย PC แต่ละตัว (proportion of variance)
- สัดส่วนข้อมูลสะสมที่อธิบายได้ เมื่อเพิ่ม PC แต่ละตัว (cumulative proportion) อย่างในตัวอย่าง จะเห็นว่า เมื่อเพิ่ม PC ตัวที่ 5 เราสามารถอธิบายข้อมูลได้ถึง 80% แล้ว แสดงว่า เราสามารถเก็บตัวแปรไว้ 5 จาก 13 ตัวแปรได้
(Note: เราควรใช้ PC เท่ากับจำนวนที่สามารถอธิบายข้อมูลได้ตั้งแต่ 80% ขึ้นไป)

วิธีที่ 2. เราสามารถดูผลลัพธ์ผ่านกราฟได้ เช่น scree plot ที่แสดงจำนวน PC และสัดส่วนข้อมูลสะสมที่อธิบายได้

ในการเริ่มสร้างกราฟ ให้เราคำนวณหาสัดส่วน variance และสัดส่วนสะสมที่อธิบายได้ก่อน:
```
# Extract variance explained
pca_var <- pca$sdev^2
pca_var_exp <- pca_var / sum(pca_var)

# Compute cumulative variance explained
cum_var_exp <- cumsum(pca_var_exp)
```
จากนั้น นำผลลัพธ์ไปสร้างกราฟ:
```
# Plot a scree plot for cumulative variance explained
plot(cum_var_exp, 
     type = "b", col = "blue", pch = 19, lwd = 2,
     main = "Cumulative Variance Explained",
     xlab = "Number of Principal Components",
     ylab = "Cumulative Variance Explained",
     ylim = c(0, 1))
```
ผลลัพธ์:

จากกราฟ เราจะได้ข้อสรุปเดียวกันกับค่าทางสถิติ นั่นคือ เมื่อเรามี PC 5 ตัว เราจะสามารถอธิบายข้อมูลได้ถึง 80% ของข้อมูลเก่า

😺 GitHub

ดู code ทั้งหมดในบทความนี้ได้ที่ GitHub

📃 References
✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

📕 ขอฝากหนังสือเล่มแรกในชีวิตด้วยนะครับ 😆

🙋 ใครที่กำลังเรียนจิตวิทยาหรือทำงานสายจิตวิทยา และเบื่อที่ต้องใช้ software ราคาแพงอย่าง SPSS และ Excel เพื่อทำข้อมูล

💪 ผมขอแนะนำ R Book for Psychologists หนังสือสอนใช้ภาษา R เพื่อการวิเคราะห์ข้อมูลทางจิตวิทยา ที่เขียนมาเพื่อนักจิตวิทยาที่ไม่เคยมีประสบการณ์เขียน code มาก่อน

ในหนังสือ เราจะปูพื้นฐานภาษา R และพาไปดูวิธีวิเคราะห์สถิติที่ใช้บ่อยกัน เช่น:
- Correlation
- t-tests
- ANOVA
- Reliability
- Factor analysis
🚀 เมื่ออ่านและทำตามตัวอย่างใน R Book for Psychologists ทุกคนจะไม่ต้องพึง SPSS และ Excel ในการทำงานอีกต่อไป และสามารถวิเคราะห์ข้อมูลด้วยตัวเองได้ด้วยความมั่นใจ

แล้วทุกคนจะแปลกใจว่า ทำไมภาษา R ง่ายขนาดนี้ 🙂‍↕️

👉 สนใจดูรายละเอียดหนังสือได้ที่ meb:

ดูรายละเอียดหนังสือ R Book for Psychologists
Share this:
X
Facebook
Like Loading…
2025-05-15

Category: Data analytics

1️⃣ Package 1. readxl

2️⃣ Package 2. XLConnect

📔 กรณีที่ 1. โหลดและสำรวจ Workbook

⬇️ กรณีที่ 2. โหลดข้อมูลจาก sheet

🖐️ กรณีที่ 3. จัดการ Sheet

➕ กรณีที่ 4. เพิ่มข้อมูลใน Sheet

💾 กรณีที่ 5. บันทึก Workbook

💪 Summary

😺 GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

Share this:

🧮 1. Comparison & Logical Operators

#️⃣ (1) Comparison Operators

♟️ (2) Logical Operators

🚦 2. Conditional Statements

🌊 3. Control Flow Statements

📦 (1) for Loop

🔁 (2) while Loop

🚄 4. Loop Control Statements

💪 5. Summary

😺 GitHub

📃 References

Share this:

💻 Install & Import

1️⃣ Basics: 5 Parametres to Know

💲 Symbols

🏦 src

🏧 auto.assign

🌲 env

📆 from, to

2️⃣ Advanced: Set Defaults

🌍 getDefaults(), setDefaults()

💹 setSymbolLookup(), getSymbolLookup()

💾 saveSymbolLookup(), loadSymbolLookup()

3️⃣ View Columns

4️⃣ Visualise Data

💪 Summary

😺 GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

Share this:

👨‍🔬 Intro to REGEX Formulas

⌨️ Regular Expression Basics

🔣 (1) Special Characters

🔣 (2) Character Classes

🔣 (3) Quantifiers

🔣 (4) Others

📃 Learn More

💪 Put It All Together

🗝️ (1) REGEXMATCH

🗝️ (2) REGEXEXTRACT

🗝️ (3) REGEXREPLACE

📚 References

Share this:

🏁 Getting Started

1️⃣ Install & Load Packages

2️⃣ Load Dataset

🔄️ Converting

1️⃣ การแปลงข้อมูลให้เป็น Datetime

2️⃣ การแปลงข้อมูลให้เป็น Time Series

🤷‍♂️ Missing Values

📈 Plotting

✂️ Subsetting

1️⃣ window()

2️⃣ Boolean Masking

3️⃣ Specific Date

🧮 Aggregating

🪟 Rolling Window

1️⃣ roll*()

2️⃣ rollapply()

➡️ Expanding Window

😎 Summary

😺 View Code in GitHub

📃 References

✅ R Book for Psychologists: หนังสือภาษา R สำหรับนักจิตวิทยา

Share this:

🚀 XGBoost คืออะไร?

💻 XGBoost ในภาษา R

**1️⃣ roll*()**