เมื่อต้องการหาความสัมพันธ์อะไรบางอย่างระหว่างตัวแปรต้นและตัวแปรตามจากข้อมูลชุดหนึ่งที่มี เรียกปัญหานี้ว่าการวิเคราะห์การถดถอย (回归, regression)

เช่นสมมุติว่ามีตัวแปร ๒ ​ตัว ตัวแปรตั้งต้น x และตัวแปรตาม z มีความสัมพันธ์แบบนี้



เนื่องจากข้อมูลนี้ตัวแปรทั้ง ๒ ดูแล้วสามารถลากเส้นตรงเพื่ออธิบายได้ ดังนั้นจึงสามารถใช้การวิเคราะห์การถดถอยเชิงเส้น (线性回归, linear regression) นั่นคือเขียน z ในรูปของ
$z = wx + b$ ..(1)

โดย w และ b เป็นค่าคงที่บางอย่างที่ต้องหา โดยอาศัยข้อมูล เช่นตัวอย่างนี้จะได้ว่า w=0.3 และ b=0.4 และจะสามารถคำนวณลากเส้นได้ผ่านกลางจุดข้อมูลแบบนี้



เรื่องของการวิเคราะห์การถดถอยเชิงเส้นได้เขียนถึงไปแล้วใน https://phyblas.hinaboshi.com/20161210

แต่ในปัญหาโดยทั่วไปแล้วยากที่จะอธิบายด้วยเส้นตรงง่ายๆแบบนั้น เช่นหากข้อมูลมีหน้าตาแบบนี้



แบบนี้แทนที่จะแค่อธิบายด้วยเส้นตรง wx+b ง่ายๆ ควรใช้ฟังก์ชันอะไรบางอย่างที่ซับซ้อนกว่านั้น

ทางหนึ่งที่สามารถทำได้ก็คือ ใช้ฟังก์ชันฐาน (基函数, basis function)

คือ z อาจเขียนอยู่ในรูปของฟังก์ชันอะไรบางอย่างของ x มาบวกกัน อาจเขียนได้ในลักษณะแบบนี้
$\begin{align} z = \sum_{i=1}^{m}w_i \phi_i(x) = w_1\phi_1(x)+w_2\phi_2(x)+\cdots+w_m\phi_m(x) \end{align}$ ..(2)

โดย ϕ_i(x) คือฟังก์ชันฐาน คือเป็นกลุ่มฟังก์ชันอะไรสักอย่างของ x ที่ต้องการนำมาใช้อธิบายค่า โดยที่ ϕ แต่ละค่านั้นจะมีคุณสมบัติตั้งฉากซึ่งกันและกัน (正交, orthogonal) คือไม่สามารถเขียนตัวนึงในรูปของตัวอื่นๆบวกกันได้

ส่วน w_i เป็นค่าน้ำหนักที่บ่งบอกว่าฟังก์ชันฐานแต่ละตัวมีความสำคัญแค่ไหน

และ m คือจำนวนฟังก์ชันฐานที่ใช้

ตัวอย่างฟังก์ชันฐานที่นิยมใช้ที่สุดก็คือ ฟังก์ชันพหุนาม นั่นคือ
$\phi_i(x) = x^{i-1}$ ..(3)

ก็จะได้ว่า
$\begin{align} z = \sum_{i=1}^{m}w_i x^{i-1} = w_1+w_2x+w_3x^2+\cdots+w_mx^{m-1} \end{align}$ ..(4)

ทีนี้เพื่อความสะดวกในการคำนวณ เราอาจเขียนใหม่ในรูปเวกเตอร์ได้เป็นแบบนี้
$\begin{align} z &= \begin{bmatrix} \phi_1(x) & \phi_2(x) & \cdots & \phi_m(x) \end{bmatrix} \begin{bmatrix} w_1 \\\\ w_2 \\\\ \vdots \\\\ w_m \end{bmatrix} \\ &= \vec{\phi}^T\vec{w} \end{align}$ ..(5)

โดย
$\begin{align} \vec{\phi}(x) = \begin{bmatrix} \phi_1(x) \\ \phi_2(x) \\ \vdots \\ \phi_m(x) \end{bmatrix} \end{align}$ ..(6)
$\begin{align} \vec{w} = \begin{bmatrix} w_1 \\ w_2 \\ \vdots \\ w_m \end{bmatrix} \end{align}$ ..(7)

สำหรับฟังก์ชันฐานพหุนาม
$\begin{align} \vec{\phi}(x) = \begin{bmatrix} 1 \\ x \\ x^2 \\ \vdots \\ x^{m-1} \end{bmatrix} \end{align}$ ..(8)

ทีนี้จะเห็นว่าค่า z ที่คำนวณได้นั้นจะขึ้นอยู่กับ w ทั้งหลาย ดังนั้นปัญหาต่อไปก่อนอื่นสิ่งที่จะต้องทำก็คือหาค่า w ที่จะทำให้คำนวณค่า z ออกมาได้ใกล้เคียงกับข้อมูลที่มีอยู่มากที่สุด

การแก้หาค่า w นั้นอาจมีอยู่หลายวิธี แต่โดยทั่วไปที่เป็นพื้นฐานง่ายที่สุดก็คือคำนวณค่าผลรวมความคลาดเคลื่อนกำลังสอง (和方差, sum of squared error, SSE) ระหว่างค่าจริงกับค่าที่คำนวณได้
$\begin{align} J(w) = \sum_{i=1}^n\left(z_{i,คำนวณ}-z_{i,จริง}\right)^2 \end{align}$ ..(9)

ทีนี้ เมื่อเราจำเป็นต้องคำนวณค่า z สำหรับค่า x ทุกตัวที่มี โดยทั่วไปแล้วทางที่สะดวกคือเขียนอยู่ในรูปของการคูณเมทริกซ์ แบบนี้
$\vec{z} = \mathbf{\Phi}\vec{w}$ ..(10)

โดย
$\begin{align} \vec{z} = \begin{bmatrix} z_1 \\ z_2 \\ \vdots \\ z_n \end{bmatrix} \end{align}$ ..(11)

โดย n คือจำนวนจุดข้อมูล

และ
$\begin{align} \mathbf{\Phi} = \begin{bmatrix} \phi_1(x_1) & \phi_2(x_1) & \cdots & \phi_m(x_1) \\ \\ \phi_1(x_2) & \phi_2(x_2) & \cdots & \phi_m(x_2) \\ \\ \vdots & \vdots & \ddots & \vdots \\ \\ \phi_1(x_n) & \phi_2(x_n) & \cdots & \phi_m(x_n) \end{bmatrix} \end{align}$ ..(12)

เมทริกซ์ของ Φ นี้มีชื่อเรียกเฉพาะว่าเมทริกซ์ออกแบบ (设计矩阵, design matrix)

แทนลงในสมการ (9) จะได้ว่า
$\begin{align} J(w) = \sum_{i=1}^n\left(\sum_{j=1}^m w_j\phi_j(x_i)-z_i\right)^2 \end{align}$ ..(13)

เป้าหมายคือต้องการหาค่า w ที่ทำให้ J ต่ำสุด โดยทั่วไปที่จุดต่ำสุดจะมีค่าความชันเป็น 0 ดังนั้นจึงหาอนุพันธ์ย่อยของ J เทียบ w และให้เท่ากับ 0 จะได้
$\begin{align} \frac{\partial J(w)}{\partial w_k} = 2\sum_{i=1}^n\phi_k(x_i)\left(\sum_{j=1}^m w_j\phi_j(x_i)-z_i\right) = 0 \end{align}$ ..(14)

แล้วก็จะได้ว่า
$\begin{align} (\mathbf{\Phi}^T\mathbf{\Phi}) \vec{w} = \mathbf{\Phi}^T \vec{z} \end{align}$ ..(15)

ซึ่งในรูปแบบนี้สามารถหาค่า w ได้ด้วยการแก้ระบบสมการ

ในไพธอนมีฟังก์ชัน np.linalg.solve ซึ่งใช้สำหรับแก้สมการแบบนี้ได้ทันทีอยู่แล้ว จึงดึงมาใช้ได้เลย เสร็จแล้วก็จะได้ค่า w ที่ต้องการออกมา

เสร็จแล้วพอได้ w มาก็นำมาหาค่า z ของข้อมูลใหม่ได้โดยสมการ (10) ต่อไป

สรุปขั้นตอนสำหรับการหาค่า z ใหม่ที่ไม่รู้คือ
1. คำนวณ Φ ของข้อมูล x ที่มี
2. ใช้ np.linalg.solve แก้หาค่า w
3. คำนวณ Φ ของข้อมูล x ใหม่ที่ต้องการหา
4. คำนวณ z จาก Φ ใหม่และ w

ลองเขียนโค้ดเพื่อใช้งานดูเลย ข้อมูลตัวอย่างที่ใช้นี้เป็นอันเดียวกับที่แสดงในภาพข้างบน ในที่นี้ใช้ x และ z เป็นตัวแปรของข้อมูลที่มี ส่วน x_ เป็นข้อมูลจุดใหม่ที่ต้องการหา z_ เป็นค่าของจุดใหม่ที่คำนวณได้ ส่วน Φ เขียนแทนด้วย phi

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(0)
n = 70 # จำนวนข้อมูล
x = np.random.uniform(0,10,n) # จุดข้อมูลที่รู้ค่า
# ค่าของแต่ละจุด ให้เป็นฟังก์ชันพหุนามที่บวกคลื่นรบกวนเข้าไป
z = ((2*x-1)**2*(x-5)**2*(x-9)**2)*np.random.normal(1,0.1,n)/1000

# ฟังก์ชันสำหรับสร้าง Φ โดยใช้ฐานเป็นพหุนาม
def than_phahunam(x,m):
    return np.array([x**i for i in range(m)]).T

m = 8 # จำนวนฐานที่จะใช้
phi = than_phahunam(x,m)
# แก้สมการหาค่า w
w = np.linalg.solve(phi.T.dot(phi),phi.T.dot(z))

x_ = np.linspace(x.min(),x.max(),201) # ข้อมูลใหม่ที่ต้องการหาค่า
phi_ = than_phahunam(x_,m)
z_ = phi_.dot(w)

plt.axes(xlabel='x',ylabel='z')
plt.scatter(x,z,c='m',edgecolor='b')
plt.plot(x_,z_,'r')
plt.show()

กรณีที่ใช้ฐานเป็นพหุนามแบบนี้อาจเรียกว่าการวิเคราะห์การถดถอยพหุนาม (多项式回归, polynomial regression)

เพื่อความสะดวกในการใช้ ลองเขียนในรูปแบบคลาสได้ดังนี้

class ThotthoiPhahunam:
    def __init__(self,m=4):
        self.m = m

    def rianru(self,x,z):
        phi = self.than(x)
        self.w = np.linalg.solve(phi.T.dot(phi),phi.T.dot(z))

    def thamnai(self,x_):
        phi_ = self.than(x_)
        return phi_.dot(self.w)

    def than(self,x):
        return x[:,None]**np.arange(self.m)

np.random.seed(0)
n = 50
x = np.random.uniform(0,20,n)
z = ((x-1)**2*(x-10)*(x-18)**2)*np.random.normal(1,0.1,n)/1000

tt = ThotthoiPhahunam(m=8)
tt.rianru(x,z)
x_ = np.linspace(x.min(),x.max(),201)
z_ = tt.thamnai(x_)

plt.axes(xlabel='x',ylabel='z')
plt.scatter(x,z,c='m',edgecolor='b')
plt.plot(x_,z_,'r')
plt.show()

ที่จริงเรื่องการวิเคราะห์การถดถอยพหุนามเคยเขียนถึงไปแล้วใน https://phyblas.hinaboshi.com/20161219

เพียงแต่ว่ารูปแบบการเขียน และวิธีการหาค่า w ก็ต่างกัน ในนั้นใช้วิธีการเคลื่อนลงตามความชัน (梯度下降法, gradient descent) แต่ในนี้ใช้วิธีการแก้สมการหาโดยตรง



ใช้ฟังก์ชันเกาส์เป็นฐาน

นอกจากพหุนามแล้ว ฐานอีกชนิดหนึ่งที่นิยมคือใช้ฟังก์ชันเกาส์ (Gauss)

ฟังก์ชันเกาส์ เป็นฟังก์ชันที่มีค่าเป็นลักษณะระฆังคว่ำ โดยมีค่าสูงสุดอยู่ที่จุดนึง เรียกว่าเป็นจุดศูนย์กลาง และยิ่งห่างจากจุดนี้ไปก็จะมีค่าต่ำลงเรื่อยๆจนเข้าใกล้ 0 นิยามดังนี้
$\begin{align} \phi(x) = \exp\left({-\gamma(x-\mu)^2}\right) \end{align}$ ..(16)

โดย μ คือจุดศูนย์กลาง ส่วน γ เป็นตัวกำหนดว่าค่าจะลงลงเร็วแค่ไหนเมื่อห่างจากใจกลาง ยิ่งค่าน้อยยิ่งลดลงช้า

เมื่อนำมาใช้เป็นฟังก์ชันฐาน ฐานเกาส์อาจนิยามตามนี้
$\begin{align} \phi_j(x_i) = \exp\left({-\gamma(x_i-x_j)^2}\right) \end{align}$ ..(17)

โดย x_j เป็นจุดข้อมูล นั่นคือใช้จุดข้อมูลแต่ละจุดเป็นศูนย์กลางของฐานแต่ละอัน

กรณีนี้เราจะมีจำนวนฐานเท่ากับจำนวนข้อมูลทั้งหมดที่ใช้

เขียนโปรแกรมได้ดังนี้

def gauss(x1,x2,gamma):
    return np.exp(-gamma*(x1-x2)**2)

n = 20
np.random.seed(4)
x = np.random.uniform(0,10,n)
z = ((x-1)**2*(x-9)**2)*np.random.normal(1,0.1,n)/1000

gamma = 1
phi = gauss(x,x[:,None],gamma)
w = np.linalg.solve(phi.T.dot(phi),phi.T.dot(z))

x_ = np.linspace(x.min(),x.max(),201)
phi_ = gauss(x,x_[:,None],gamma)
z_ = phi_.dot(w)

plt.axes(xlabel='x',ylabel='z')
plt.scatter(x,z,c='m',edgecolor='b')
plt.plot(x_,z_,'r')
plt.show()

ผลที่ได้จะพบว่าออกมาแปลกๆ ที่เป็นแบบนี้เพราะว่าเมื่อจำนวนฐานเท่ากับจำนวนข้อมูลจะทำให้ฟังก์ชันที่ได้สามารถลากผ่านจุดข้อมูลทั้งหมด

ซึ่งที่จริงแล้วข้อมูลโดยทั่วไปจะมีความไม่แน่นอนปนอยู่ ดังนั้นการฝืนให้เส้นลากผ่านจุดทั้งหมดนั้นจึงไม่ใช่เรื่องที่ดี มีแต่จะทำให้เกิดการเรียนรู้เกิน (过学习, overlearning) คือการที่ผลการคำนวณปรับตัวเข้ากับข้อมูลที่มีอยู่มากเกินไปแต่ไม่สามารถใช้อธิบายทุกอย่างตามธรรมชาติจริงๆได้

เมื่อเป็นแบบนี้จึงจำเป็นต้องมีการเรกูลาไรซ์ (正规化, regularize) คือการทำให้การคำนวณไม่ยึดติดกับข้อมูลที่มีมากเกินไป

เกี่ยวกับแนวคิดของการเรกูลาไรซ์เคยเขียนถึงไปแล้วใน https://phyblas.hinaboshi.com/20170928

ในที่นี้ใช้เรกูลาไรซ์แบบ l2 สมการ (13) ก็จะแก้เป็น
$\begin{align} J(w) = \sum_{i=1}^n\left(\sum_{j=1}^m w_j\phi_j(x_i)-z_i\right)^2 + \lambda\sum_{j=1}^{m}w_j^2 \end{align}$ ..(18)

ที่เปลี่ยนไปคือมีเพิ่มพจน์ทางขวาบวกเข้ามา โดย λ คือขนาดของการเรกูลาไรซ์

จากนั้นก็ให้ความชันเท่ากับ 0 เช่นเคย
$\begin{align} \frac{\partial J(w)}{\partial w_k} = 2\sum_{i=1}^n\phi_k(x_i)\left(\sum_{j=1}^m w_j\phi_j(x_i)-z_i\right) + 2\lambda w_k = 0 \end{align}$ ..(19)

แล้วก็จะได้
$\begin{align} (\mathbf{\Phi}^T\mathbf{\Phi}+\lambda I) \vec{w} = \mathbf{\Phi}^T \vec{z} \end{align}$ ..(20)

โดย I คือเมทริกซ์เอกลักษณ์ ค่า λ ที่เพิ่มเข้าไปจะมีผลทำให้ค่าแนวทแยงของเมทริกซ์มากขึ้น

สามารถเขียนโค้ดแก้ตรงที่คำนวณ w ใหม่ เป็นแบบนี้

l = 0.1
w = np.linalg.solve(phi.T.dot(phi) + np.eye(len(x))*l,phi.T.dot(z))

แล้วผลที่ได้จะกลายเป็นแบบนี้



จะเห็นว่าได้เส้นกราฟออกมาเรียบและดูสมเหตุสมผลขึ้น



สุดท้ายนี้ เพื่อความสะดวกในการใช้งานลองมาสร้างเป็นคลาสขึ้นมาดู

class ThotthoiThanGauss:
    def __init__(self,gamma=1,l=0.1):
        self.gamma = gamma
        self.l = l

    def rianru(self,x,z):
        phi = self.gauss(x,x[:,None])
        self.w = np.linalg.solve(phi.T.dot(phi)+self.l*np.eye(len(x)),phi.T.dot(z))
        self.x = x

    def thamnai(self,x_):
        phi_ = self.gauss(self.x,x_[:,None])
        return phi_.dot(self.w)

    def gauss(self,x1,x2):
        return np.exp(-self.gamma*(x1-x2)**2)

ข้อควรสังเกตอย่างหนึ่งคือ กรณีนี้จำเป็นต้องบันทึกตำแหน่งของ x ที่ป้อนเข้ามาตอนแรกไว้ด้วย เพราะต้องใช้ตำแหน่งใจกลางของฐานเกาส์ในการคำนวณ

ขอยกตัวอย่างการใช้ด้วยการเปรียบเทียบความแตกต่างระหว่างกรณีที่ λ มีค่าต่างกัน

n = 30
np.random.seed(5)
x = np.random.uniform(0,30,n)
z = ((x-5)**2*(x-15)*(x-28))*np.random.normal(1,0.2,n)/1000
x_ = np.linspace(x.min(),x.max(),301)
ll = [0,0.1,1,10]
plt.figure(figsize=[6.5,5])
for i in range(4):
    l = ll[i]
    tt = ThotthoiThanGauss(gamma=0.1,l=l)
    tt.rianru(x,z)
    z_ = tt.thamnai(x_)
    plt.subplot(2,2,i+1)
    plt.scatter(x,z,c='m',edgecolor='b')
    plt.title('$\\lambda$=%.1f'%l)
    plt.plot(x_,z_,'r')
plt.tight_layout()
plt.show()

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



ทั้งหมดนี้เป็นพื้นฐานของการประยุกต์ใช้ฟังก์ชันฐานในปัญหาการวิเคราะห์การถดถอย

การใช้ฟังก์ชันฐานยังเป็นพื้นฐานสำคัญที่จะต่อยอดไปยังเรื่องอื่นในเทคทิคการเรียนรู้ของเครื่อง เช่นการถดถอยเชิงเส้นแบบเบส์ (贝叶斯线性回归, Bayesian Linear Regression) รวมทั้งการใช้ลูกเล่นเคอร์เนล (核技巧, kernel trick)

เกี่ยวกับลูกเล่นเคอร์เนลอ่านได้ใน https://phyblas.hinaboshi.com/20180724



อ้างอิง

~ เกี่ยวกับเรา ~

สารบัญ

รวมคำแปลวลีเด็ดจากญี่ปุ่น

ภาษา python มอดูลต่างๆ -- numpy -- matplotlib -- pandas -- manim -- opencv -- pyqt -- pytorch การเรียนรู้ของเครื่อง -- โครงข่าย      ประสาทเทียม

ภาษา javascript

ภาษา mongol

ภาษาศาสตร์

maya

ความน่าจะเป็น

บันทึกในญี่ปุ่น

บันทึกในจีน -- บันทึกในปักกิ่ง -- บันทึกในฮ่องกง -- บันทึกในมาเก๊า

บันทึกในไต้หวัน

บันทึกในยุโรปเหนือ

บันทึกในประเทศอื่นๆ

qiita

บทความอื่นๆ

บทความแบ่งตามหมวด

==เลือกหมวด== ดาราศาสตร์ คณิตศาสตร์ -ความน่าจะเป็น คอมพิวเตอร์ -เขียนโปรแกรม --python ---numpy ---scipy ---matplotlib ---pandas ---manim ---pyqt ---sklearn ---pytorch ---mayapython --ruby --javascript --dart --MATLAB --SQL --regex --opencv -shell -3D --maya --MMD -microsoft_office -pdf -ปัญญาประดิษฐ์ --โครงข่ายประสาทเทียม --สเตเบิลดิฟฟิวชัน ---comfyui -การสุ่ม ภาษาศาสตร์ -ตัวอักษร -เรียนภาษา -หลักเกณฑ์การทับศัพท์ -ภาษาจีน --ภาษาจีนกลาง -ภาษาญี่ปุ่น -ภาษามองโกล -ภาษาลาว -ภาษาเขมร ประวัติศาสตร์ -ประวัติศาสตร์จีน -ประวัติศาสตร์ญี่ปุ่น ปรัชญา ประเทศจีน -จีนแผ่นดินใหญ่ --ปักกิ่ง --เทียนจิน --เหลียวหนิง --เหอเป่ย์ --เหอหนาน --ซานตง --ซานซี --อานฮุย --เจ้อเจียง --หูเป่ย์ --หูหนาน --ฝูเจี้ยน --กวางตุ้ง ---แต้จิ๋ว --ยูนนาน --ซินเจียง -ฮ่องกง -มาเก๊า -ไต้หวัน --ไทเป --จีหลง --เถาหยวน --ซินจู๋ --เหมียวลี่ --ไถจง --จางฮว่า --หยวินหลิน --เจียอี้ --ไถหนาน --เกาสยง --ผิงตง --อี๋หลาน ประเทศญี่ปุ่น -ฮกไกโด -อาโอโมริ -อิวาเตะ -มิยางิ -อากิตะ -ยามางาตะ -ฟุกุชิมะ -อิบารากิ -โทจิงิ -กุมมะ -ไซตามะ -จิบะ -โตเกียว -คานางาวะ -นีงาตะ -โทยามะ -อิชิกาวะ -ฟุกุอิ -ยามานาชิ -นางาโนะ -กิฟุ -ชิซึโอกะ -ไอจิ -มิเอะ -ชิงะ -เกียวโต -โอซากะ -เฮียวโงะ -นาระ -วากายามะ -โอกายามะ -ฮิโรชิมะ -ยามางุจิ -ฟุกุโอกะ -ซางะ -นางาซากิ -คุมาโมโตะ -โออิตะ ต่างแดน -อุษาคเนย์ --กัมพูชา --พม่า --สิงคโปร์ -ยุโรป --สวีเดน --เดนมาร์ก --ฟินแลนด์ -ออสเตรเลีย ท่องเที่ยว -มหาวิทยาลัย -พิพิธภัณฑ์ --พิพิธภัณฑ์วิทยาศาสตร์ --หอศิลป์ -สวนสัตว์ --พิพิธภัณฑ์สัตว์น้ำ -ท้องฟ้าจำลอง -ตึกระฟ้า -ปราสาท☑ --ปราสาทญี่ปุ่น --ปราสาทขอม --ปราสาทยุโรป -ศาสนสถาน --วัด --ศาลเจ้า --โบสถ์ --มัสยิด -สุสาน -มรดกโลก -ทะเล -ทะเลสาบ -ภูเขา -หิมะ -ดอกซากุระ -แมว -รถไฟ -เรือ -ตลาดกลางคืน -งานเทศกาล -ที่ระลึกภัยพิบัติ -ตามรอย บันเทิง -เกม --อาเตอลีเย --โปเกมอน --caligula --vn -อนิเมะ -มังงะ -นิยาย -เพลง --เพลงอนิเมะ --เพลงเกม เรื่องแต่ง บันทึก

　 ค้นหาบทความ

　 บทความล่าสุด

　 บทความแนะนำ

บทความแต่ละเดือน

ค้นบทความเก่ากว่านั้น

ไทย

日本語

中文