มาร์คอฟเชนมอนติคาร์โล (MCMC): สุ่มจากการกระจายความน่าจะเป็นที่ซับซ้อน

การกระจายความน่าจะเป็นบางแบบซับซ้อนเกินกว่าจะทำงานด้วยโดยตรง

เมื่อต้องสร้างแบบจำลองข้อมูลโลกจริง คณิตศาสตร์มักไปต่อไม่ไหวก่อนจะได้ข้อสรุปที่ใช้การได้ หลายครั้งอินทิกรัลดูเหมือนจัดการได้บนกระดาษ แต่กลับแก้ไม่ออกทันทีที่เพิ่มตัวแปรแฝงสักสองสามตัว เรื่องนี้พบได้บ่อยในอนุมานแบบเบย์เชียน ซึ่งการกระจายภายหลังจะผสานความเชื่อก่อนหน้าเข้ากับข้อมูลที่สังเกตได้ และผลลัพธ์ที่ได้มักไม่สามารถสรุปด้วยสูตรง่าย ๆ

แนวคิดพื้นฐานของ Markov Chain Monte Carlo คือ แทนที่จะลงมือกับคณิตศาสตร์โดยตรง MCMC จะสำรวจการกระจายผ่านการจำลอง โดยดึงตัวอย่างที่สะท้อนรูปร่างของมันได้ โดยไม่ต้องคำนวณแบบเต็มรูป

บทความนี้จะอธิบายแนวคิดแกนหลักของ MCMC ไล่ดูอัลกอริทึมที่พบบ่อยที่สุด และสาธิตวิธีใช้งานใน Python

ต้องการทบทวนคณิตศาสตร์สำหรับ Python หรือไม่ อ่านบล็อกโพสต์ของเราเรื่อง Demystifying Mathematical Concepts for Deep Learning เพื่อดูตัวอย่างคณิตศาสตร์ใน Numpy

Markov Chain Monte Carlo คืออะไร

Markov Chain Monte Carlo (MCMC) คือกลุ่มอัลกอริทึมที่สร้างตัวอย่างจากการกระจายความน่าจะเป็น แม้เมื่อการกระจายนั้นซับซ้อนเกินกว่าจะทำงานด้วยโดยตรง

ชื่อเรียกแบ่งได้เป็นสองส่วน ส่วนที่เป็นมาร์คอฟเชนควบคุมวิธีที่อัลกอริทึมเคลื่อนผ่านสถานะที่เป็นไปได้ แต่ละก้าวขึ้นกับสถานะปัจจุบันเท่านั้น ไม่ขึ้นกับเส้นทางในอดีต ส่วนมอนติคาร์โลหมายถึงการใช้การสุ่มตัวอย่างเพื่อประมาณค่าที่สนใจ

เมื่อรวมกัน MCMC จะสร้างสายโซ่ของตัวอย่างสุ่มซึ่งเมื่อเวลาผ่านไปจะสะท้อนรูปร่างของการกระจายเป้าหมาย เป็นเทคนิคการสุ่มตัวอย่างก่อนอื่น ไม่ได้แก้คณิตศาสตร์อย่างแม่นยำ แต่ประมาณด้วยการจำลอง

ทำไมต้องใช้ Markov Chain Monte Carlo

ปัญหาของการกระจายข้อมูลโลกจริงคือมันไม่ได้สวยงามเหมือนในตำรา

ในอนุมานแบบเบย์เชียน มักต้องคำนวณการกระจายภายหลัง ซึ่งเป็นความน่าจะเป็นที่อัปเดตของพารามิเตอร์หลังเห็นข้อมูล สูตรบนกระดาษดูง่าย: คูณ prior กับ likelihood แล้วหารด้วย marginal likelihood พจน์สุดท้ายนั้นต้องอินทิเกรตเหนือค่าพารามิเตอร์ทั้งหมด ในมิติสูง อินทิกรัลนี้แทบคำนวณไม่ได้จริง

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

• แบบจำลองเชิงลำดับชั้นแบบเบย์เชียน: ที่พารามิเตอร์ขึ้นกับพารามิเตอร์อื่น ทำให้เกิดการกระจายซ้อนชั้นที่อินทิเกรตไม่ได้แบบสะอาด
• การกระจายภายหลังมิติสูง: จำนวนพารามิเตอร์มากเกินกว่าจะประมาณแบบ grid ได้อย่างมีประสิทธิภาพ
• priors ที่ไม่เป็นคอนจูเกต: ที่ prior และ likelihood ไม่ผสานกันเป็นการกระจายที่รู้จัก

MCMC เป็นทางออกที่ดีในสถานการณ์เหล่านี้ แทนที่จะคำนวณการกระจาย ก็สุ่มตัวอย่างจากมัน ตัวอย่างเหล่านี้ให้ทุกอย่างที่ต้องการโดยไม่ต้องแก้อินทิกรัล

สองแนวคิดเบื้องหลัง MCMC

MCMC ผสานสองแนวคิดที่ต่างก็เรียบง่ายเมื่อแยกกัน แต่ทรงพลังเมื่อรวมกัน มาดูกันทีละส่วน

มาร์คอฟเชน

มาร์คอฟเชนคือชุดของสถานะที่แต่ละก้าวขึ้นกับตำแหน่งปัจจุบันเท่านั้น

อดีตที่ผ่านมาไม่สำคัญ มีเพียงสถานะปัจจุบันที่กำหนดว่าจะไปไหนต่อ คุณสมบัติ “ไร้ความทรงจำ” นี้ ซึ่งทางการเรียกว่าสมบัติมาร์คอฟ ทำให้คณิตศาสตร์จัดการได้และอัลกอริทึมใช้งานได้จริง

โซ่จะเคลื่อนผ่านสถานะทีละก้าว และด้วยการตั้งค่าที่เหมาะสม ในที่สุดมันจะนิ่งที่การกระจายสถานีแนรี ซึ่งเป็นแบบแผนคงที่ที่ความน่าจะเป็นของแต่ละสถานะไม่เปลี่ยนแปลง การกระจายที่เสถียรนี้เองคือสิ่งที่ MCMC ออกแบบมาให้ทำงานด้วย

วิธีมอนติคาร์โล

วิธีมอนติคาร์โลใช้การสุ่มตัวอย่างเพื่อประมาณค่าที่คำนวณโดยตรงได้ยาก

แนวคิดคือสุ่มตัวอย่างจากการกระจายให้มากพอ แล้วประมาณค่าเฉลี่ย ความแปรปรวน หรือคุณสมบัติอื่น ๆ จากตัวอย่างเหล่านั้น ยิ่งมีตัวอย่างมาก ผลประมาณก็ยิ่งใกล้ค่าจริง

โดยลำพัง วิธีมอนติคาร์โลต้องสุ่มโดยตรงจากการกระจาย ซึ่งเป็นปัญหาที่เราพยายามแก้ มาร์คอฟเชนเข้ามาช่วยส่วนนี้

MCMC ทำงานอย่างไร

MCMC เป็นลูปที่มีการตัดสินใจง่าย ๆ ในแต่ละก้าว

1. เริ่มจากสถานะเริ่มต้น: เลือกจุดเริ่มในพื้นที่พารามิเตอร์ แม้จะเป็นการเดาคร่าว ๆ ก็ตาม
2. เสนอสถานะใหม่: ใช้กลไกการเสนอเพื่อชี้ว่าควรขยับไปที่ใดต่อไป
3. ตัดสินใจว่าจะยอมรับหรือไม่: เปรียบเทียบสถานะที่เสนอและสถานะปัจจุบัน แล้วรับหรือปฏิเสธตามความน่าจะเป็นภายใต้การกระจายเป้าหมาย
4. ทำซ้ำ: รันลูปนี้หลายพันครั้ง เก็บตัวอย่างระหว่างทาง
5. ใช้ตัวอย่างที่ได้: ประมาณค่าเฉลี่ย ช่วงความน่าเชื่อถือ หรือปริมาณอื่น ๆ ที่สนใจจากตัวอย่างที่เก็บมา

ขั้นตอนรับ/ปฏิเสธคือจุดที่เกิด “เวทมนตร์”

ด้วยการยอมรับสถานะที่ดีกว่าบ่อยกว่าสถานะที่แย่กว่า โซ่จะโน้มไปยังบริเวณที่มีความน่าจะเป็นสูง โดยไม่ต้องคำนวณการกระจายแบบเต็มรูป

ตัวอย่างช่วงต้นขึ้นกับจุดเริ่มต้น จึงมักถูกทิ้ง หลังจากวนซ้ำมากพอ โซ่จะลืมจุดเริ่ม และตัวอย่างที่เหลือจะสะท้อนรูปร่างจริงของการกระจายเป้าหมาย

การกระจายเป้าหมายและการกระจายสถานีแนรี

MCMC ถูกสร้างมาเพื่อสร้างตัวอย่างจากการกระจายเป้าหมายที่ไม่สามารถสุ่มตรง ๆ ได้

การกระจายเป้าหมายคือสิ่งที่ต้องการเรียนรู้ โดยมากคือการกระจายภายหลังในเบย์เชียน ทราบรูปร่างได้ถึงค่าคงที่ทำให้เป็นปกติ แต่คำนวณค่าคงที่นั้นตรง ๆ ไม่ได้ ซึ่ง MCMC ก็ไม่ต้องใช้มัน

อัลกอริทึม MCMC ทุกตัวออกแบบให้โซ่มาร์คอฟของมันมีการกระจายสถานีแนรีที่ตรงกับการกระจายเป้าหมาย การกระจายสถานีแนรีคือการกระจายที่โซ่จะนิ่งลงหลังจากเดินไปมากพอ

ปล่อยให้โซ่รันต่อไป มันจะเริ่มผลิตตัวอย่างที่ดูเหมือนการสุ่มจากการกระจายเป้าหมายของคุณ โดยเลี่ยงการอินทิเกรตไปได้

อัลกอริทึม MCMC ยอดนิยม

มีอัลกอริทึม MCMC อยู่หลายแบบที่พบในการใช้งานจริง ทั้งหมดทำตามลูปแกนเดียวกัน แต่ต่างกันที่วิธีเสนอสถานะใหม่และการใช้ข้อมูลเกี่ยวกับการกระจายเป้าหมาย

อัลกอริทึมเมโทรโปลิส

เมโทรโปลิสเป็นอัลกอริทึม MCMC ที่ง่ายที่สุดและเป็นจุดเริ่มต้นของแนวทางนี้

ในแต่ละก้าว จะเสนอสถานะใหม่ด้วยการเติมสัญญาณรบกวนแบบสุ่มลงในสถานะปัจจุบัน หากสถานะที่เสนอมีความน่าจะเป็นสูงกว่า ภายใต้การกระจายเป้าหมาย จะถูกยอมรับเสมอ หากต่ำกว่า จะยอมรับด้วยความน่าจะเป็นแปรผันตามอัตราส่วนของสองความน่าจะเป็น มิฉะนั้นโซ่จะอยู่ที่เดิม

กลไกรับ/ปฏิเสธนี้ทำให้โซ่อยู่ในบริเวณความน่าจะเป็นสูงนานขึ้น โดยไม่ต้องคำนวณการกระจายเต็มรูป

เมโทรโปลิสใช้การกระจายข้อเสนอแบบสมมาตร หมายความว่ามีโอกาสเสนอไปทุกทิศเท่า ๆ กัน ซึ่งมักเริ่มมีปัญหาเมื่อแบบจำลองใหญ่ขึ้น

อัลกอริทึมเมโทรโปลิส–เฮสติงส์

เมโทรโปลิส–เฮสติงส์ (MH) ทั่วไปกว่าเมโทรโปลิส โดยอนุญาตให้ใช้การกระจายข้อเสนอแบบไม่สมมาตร

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

วิธี MCMC สมัยใหม่ส่วนใหญ่เป็นส่วนขยายของ MH หรือสร้างบนหลักการเดียวกัน ดังนั้นหากเข้าใจเมโทรโปลิส–เฮสติงส์ ก็เข้าใจรากฐานของสาขานี้

กิบส์แซมปลิง

กิบส์แซมปลิงอัปเดตทีละตัวแปร แทนที่จะเสนอสถานะใหม่สำหรับทุกพารามิเตอร์พร้อมกัน

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

วิธีนี้เลี่ยงขั้นตอนรับ/ปฏิเสธโดยสิ้นเชิง เพราะการสุ่มแบบมีเงื่อนไขยอมรับเสมอ เหมาะเมื่อการกระจายร่วมทั้งหมดสุ่มยากแต่การมีเงื่อนไขทำได้ ซึ่งพบได้บ่อยในแบบจำลองเชิงลำดับชั้นแบบเบย์เชียน

แฮมิลโทเนียนมอนติคาร์โล

แฮมิลโทเนียนมอนติคาร์โล (HMC) เป็นอัลกอริทึมแรกที่ทำให้อนุมานแบบเบย์เชียนสมัยใหม่ใช้ได้จริงในสเกลใหญ่

แทนที่จะเสนอแบบสุ่ม HMC ใช้ข้อมูลกราเดียนต์จากการกระจายเป้าหมายเพื่อเสนอสถานะที่ไกลจากตำแหน่งปัจจุบันแต่ยังมีโอกาสถูกยอมรับสูง เคลื่อนผ่านพื้นที่พารามิเตอร์ได้ดีกว่าวิธีเดินแบบสุ่ม มีการปฏิเสธน้อยลงและสำรวจการกระจายมิติสูงได้ดีกว่า

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

HMC เป็นกลไกเบื้องหลัง Stan ซึ่งเป็นแพลตฟอร์มโปรแกรมมิงเชิงความน่าจะเป็นที่ใช้กันอย่างแพร่หลาย NUTS (No-U-Turn Sampler) ซึ่งเป็นส่วนขยายแบบปรับตัวของ HMC ที่ใช้ใน PyMC ช่วยตัดความจำเป็นในการปรับขนาดก้าวและจำนวนก้าวด้วยตนเอง

MCMC ในสถิติเบย์เชียน

หากมีพื้นที่หนึ่งที่ MCMC ส่งผลมากที่สุด ก็คืออนุมานแบบเบย์เชียน

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

MCMC ตัดขั้นตอนนี้ออกโดยสิ้นเชิง เพียงประเมินค่าการกระจายภายหลังที่ยังไม่ทำให้เป็นปกติ ณ จุดใดจุดหนึ่ง แล้วปล่อยให้โซ่ทำงานต่อ

ตัวอย่างง่าย ๆ สมมติจะประมาณอคติของเหรียญ เริ่มจากความเชื่อก่อนหน้าว่าเหรียญน่าจะยุติธรรม จากนั้นสังเกตผลการโยน สำหรับแบบจำลองเหรียญง่าย ๆ การกระจายภายหลังมีรูปแบบปิด แต่หากเพิ่มโครงสร้างเชิงลำดับชั้น เช่น ประมาณอคติของเหรียญร้อยเหรียญพร้อมกัน จะคำนวณตรง ๆ ไม่ได้

ด้วย MCMC เพียงรันโซ่ เก็บตัวอย่างจากการกระจายภายหลัง แล้วใช้ตัวอย่างเหล่านั้นคำนวณสิ่งที่ต้องการ

ทำความเข้าใจ Burn-In การบรรจบ และการมิกซ์

สามแนวคิดนี้มักทำให้นักวิทยาศาสตร์ข้อมูลที่ใหม่กับ MCMC สับสน หากเข้าใจผิด จะยังได้ผลลัพธ์ แต่ไม่รู้ว่าทำไมจึงไม่น่าเชื่อถือ

Burn-in

เมื่อโซ่มาร์คอฟเริ่มทำงาน มันยังไม่รู้ว่าบริเวณความน่าจะเป็นสูงของการกระจายเป้าหมายอยู่ที่ไหน

ตัวอย่างช่วงต้นได้รับอิทธิพลจากจุดเริ่มต้น ไม่ใช่จากการกระจายเป้าหมาย Burn-in คือการทิ้งตัวอย่างเหล่านี้ โดยรันโซ่ไปก่อนจำนวนหนึ่ง ทิ้งตัวอย่างนั้น แล้วเก็บเฉพาะสิ่งที่ตามมาหลังโซ่มีเวลาหาตำแหน่งที่ดี

ไม่มีเกณฑ์ตายตัวว่าควร burn-in นานเท่าไร ขึ้นกับแบบจำลอง จุดเริ่มต้น และคุณภาพการมิกซ์ ในทางปฏิบัติ มักวินิจฉัยด้วยภาพ trace plot มากกว่ากำหนดตัวเลขล่วงหน้า

การบรรจบ

การบรรจบหมายถึงโซ่ไม่ถูกอิทธิพลจากจุดเริ่มอีกต่อไป และกำลังดึงตัวอย่างที่สะท้อนการกระจายเป้าหมายแล้ว

โซ่ที่ยังไม่บรรจบจะให้ตัวอย่างมีอคติ ค่าเฉลี่ยที่คำนวณได้จะไม่ตรงกับค่าเฉลี่ยภายหลังที่แท้จริง แต่สะท้อนบริเวณที่โซ่ติดค้างอยู่แทน

การบรรจบเป็นสิ่งที่ประเมินภายหลังด้วยวิธีวินิจฉัย การรันหลายโซ่จากจุดเริ่มต่างกันและตรวจว่าผลสอดคล้องกันหรือไม่ เป็นวิธีที่เชื่อถือได้ในการจับความล้มเหลวของการบรรจบ

การมิกซ์

โซ่ที่บรรจบแล้วแต่มิกซ์ไม่ดีก็ยังเป็นปัญหา

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

การมิกซ์ที่แย่มักเห็นได้ใน trace plot ที่ดูเหมือนแม่น้ำไหลเอื่อย แทนที่จะเป็นแถบแนวนอนมีนอยส์ เมื่อเห็นเช่นนั้น ตัวสุ่มต้องปรับแต่ง ใช้การกระจายข้อเสนอที่ดีกว่า หรือเปลี่ยนอัลกอริทึมไปเลย

กราฟเปรียบเทียบการมิกซ์

ประเมินผลลัพธ์ของ MCMC อย่างไร

ต่อไปนี้คือสี่วิธีในการประเมิน MCMC และคำอธิบายว่าใช้เมื่อใด

Trace plots

Trace plot แสดงค่าที่สุ่มได้ของพารามิเตอร์ในแต่ละรอบ เป็นสิ่งแรกที่ควรดูหลังรัน MCMC

Trace plot ที่สุขภาพดีจะดูเหมือนสัญญาณรบกวนขาวรอบค่าเฉลี่ยคงที่ ไม่ควรเห็นแนวโน้ม ช่วงแบนยาว หรือการเลื่อนไปช้า ๆ หากเห็นโซ่เตร็ดเตร่หรือค้างนานในบริเวณเดียว นั่นคือปัญหาการมิกซ์ และตัวอย่างไม่น่าเชื่อถือ

ภาพ Trace plots

อัตสหสัมพันธ์

ตัวอย่างจาก MCMC ไม่ได้เป็นอิสระต่อกันทั้งหมด แต่ละตัวอย่างได้รับอิทธิพลจากตัวอย่างก่อนหน้า อัตสหสัมพันธ์วัดความแข็งแรงของสหสัมพันธ์ข้ามรอบ

อัตสหสัมพันธ์สูงหมายความว่าตัวอย่างให้ข้อมูลน้อยกว่าที่จำนวนระบุ ตัวอย่าง 2,000 ตัวอย่างที่มีสหสัมพันธ์อาจให้ข้อมูลพอ ๆ กับ 200 ตัวอย่างที่เป็นอิสระ ไลบรารี MCMC ส่วนใหญ่มีกราฟอัตสหสัมพันธ์เพื่อดูว่าความสัมพันธ์ลดลงเร็วเพียงใดเมื่อระยะห่างของตัวอย่างเพิ่มขึ้น

ภาพกราฟอัตสหสัมพันธ์

ขนาดตัวอย่างเชิงประสิทธิผล

ขนาดตัวอย่างเชิงประสิทธิผล (ESS) แปลงอัตสหสัมพันธ์เป็นตัวเลขที่ใช้จริง: โซ่ของคุณเทียบเท่ากับตัวอย่างอิสระกี่ตัว

หากดึง 5,000 ตัวอย่างแต่ ESS เท่ากับ 200 คุณกำลังทำงานด้วยพลังเชิงสถิติเท่ากับการจับฉลากอิสระ 200 ตัวอย่าง ESS ต่ำหมายถึงต้องรันโซ่นานขึ้น ปรับจูนตัวสุ่ม หรือทั้งสองอย่าง โดยทั่วไปมุ่งที่ ESS อย่างน้อยสองสามร้อยต่อพารามิเตอร์ก่อนเชื่อถือค่าประมาณ

ตัวชี้วัดการบรรจบ

เมื่อรันหลายโซ่ สามารถทดสอบอย่างเป็นทางการว่าโซ่บรรจบสู่การกระจายเดียวกันหรือไม่ ตัววินิจฉัยเกลแมน–รูบิน รายงานเป็น R-hat เปรียบเทียบความแปรปรวนภายในโซ่กับระหว่างโซ่

R-hat ใกล้ 1.0 หมายถึงโซ่เห็นพ้อง ซึ่งเป็นสัญญาณที่ดี ค่ามากกว่า 1.01 หรือ 1.05 (ขึ้นกับเกณฑ์ของไลบรารี) บ่งชี้ว่าโซ่ยังไม่บรรจบและต้องรันนานขึ้น ไลบรารีสมัยใหม่อย่าง PyMC คำนวณ R-hat อัตโนมัติและเตือนเมื่อค่าสูงเกินไป

MCMC ใน Python

Python มีไลบรารีสำหรับ MCMC อยู่สองสามตัว โดยมีแนวคิดต่างกัน

• PyMC: ตัวเลือกที่เข้าถึงง่ายที่สุดสำหรับนักวิทยาศาสตร์ข้อมูลส่วนใหญ่ กำหนดแบบจำลองด้วยไวยากรณ์ Python และ PyMC จัดการการสุ่มเบื้องหลังด้วย NUTS (เวอร์ชันปรับตัวของ HMC) เป็นเครื่องมือหลักสำหรับแบบจำลองแบบเบย์เชียนใน Python
• CmdStanPy / PyStan: อินเทอร์เฟซ Python สู่ Stan ซึ่งเป็นภาษาสำหรับโปรแกรมมิงเชิงความน่าจะเป็นโดยเฉพาะ Stan คอมไพล์แบบจำลองเป็น C++ ก่อนรัน จึงเร็วและเหมาะกับแบบจำลองเชิงลำดับชั้นที่ซับซ้อน แต่แลกมาด้วยเส้นโค้งการเรียนรู้ที่ชันกว่า
• NumPy: ไม่มีไลบรารี MCMC โดยเฉพาะ แต่สามารถเขียนอัลกอริทึมอย่างเมโทรโปลิส–เฮสติงส์เองได้ มีประโยชน์เพื่อเข้าใจสิ่งที่เกิดขึ้นเบื้องหลังก่อนขยับไปใช้เครื่องมือระดับสูง

สำหรับงานส่วนใหญ่ เริ่มที่ PyMC จะเหมาะ บทความนี้จะใช้เครื่องมือนี้ ดังนั้นหากจะทำตาม โปรดติดตั้งไลบรารีก่อน:

pip install pymc

เพื่อให้ง่าย จะยึดตัวอย่างพื้นฐาน คือการประมาณอคติของเหรียญจากลำดับการโยน

ขั้นที่ 1: กำหนดแบบจำลอง

import pymc as pm
import numpy as np

# 1 = heads, 0 = tails
observed_flips = np.array([1, 0, 1, 1, 0, 1, 1, 1, 0, 1])

with pm.Model() as coin_model:
    # Prior: we believe the coin is probably fair
    bias = pm.Beta("bias", alpha=2, beta=2)

    # Likelihood: observed flips given the bias
    flips = pm.Bernoulli("flips", p=bias, observed=observed_flips)

ค่า prior pm.Beta เข้ารหัสความเชื่ออ่อน ๆ ว่าเหรียญยุติธรรม ส่วน likelihood pm.Bernoulli เชื่อมแบบจำลองเข้ากับข้อมูลที่สังเกตได้

ขั้นที่ 2: รันตัวสุ่ม

with coin_model:
    trace = pm.sample(2000, tune=1000, return_inferencedata=True)

ผลลัพธ์จากการรันตัวสุ่ม

tune กำหนดจำนวนก้าว burn-in ซึ่งจะถูกทิ้ง sample ดึงตัวอย่างภายหลัง 2000 ตัวต่อโซ่หลังการจูน

ขั้นที่ 3: ตรวจสอบการกระจายภายหลัง

import arviz as az

az.plot_trace(trace, var_names=["bias"])
az.summary(trace, var_names=["bias"])

Trace plot ของแบบจำลองและผลสรุป

az.summary() ให้ค่าเฉลี่ยภายหลัง ส่วนเบี่ยงเบนมาตรฐาน และ R-hat ของแต่ละพารามิเตอร์ หาก R-hat ใกล้ 1.0 แสดงว่าโซ่บรรจบ az.plot_trace() วาด trace และการกระจายภายหลังเคียงกันสำหรับแต่ละพารามิเตอร์

สำหรับชุดข้อมูลนี้ — หัว 7 จาก 10 ครั้ง — ค่าเฉลี่ยภายหลังเท่ากับ 0.642 ด้วยส่วนเบี่ยงเบนมาตรฐาน 0.124 สะท้อนหลักฐานจากข้อมูลขณะยังคงใกล้ prior เหรียญยุติธรรม R-hat เท่ากับ 1.00 และ ESS มากกว่า 2000 มาก แสดงว่าโซ่บรรจบและตัวอย่างเชื่อถือได้

ข้อผิดพลาดที่พบบ่อยกับ MCMC

MCMC รันได้ง่าย แต่ก็ใช้ผิดได้ง่ายเช่นกัน ต่อไปนี้คือข้อผิดพลาดที่พบบ่อย

• สรุปว่าโซ่บรรจบเร็วเกินไป: แค่ตัวสุ่มรันจบโดยไม่เออเรอร์ไม่ได้แปลว่าบรรจบ ตรวจสอบ R-hat และ trace plot เสมอก่อนเชื่อผล โซ่อาจดูนิ่งแต่ยังติดอยู่ในบริเวณหนึ่งของพื้นที่พารามิเตอร์
• เมิน burn-in: ตัวอย่างช่วงแรกเอนเอียงไปทางจุดเริ่ม ไม่ใช่การกระจายเป้าหมาย ทิ้งมันไป ไลบรารีส่วนใหญ่จัดการผ่านพารามิเตอร์ tune อยู่แล้ว แต่ควรตรวจว่าไม่ได้เผลอรวมตัวอย่างเหล่านั้นในงานวิเคราะห์
• การมิกซ์แย่: โซ่ที่เคลื่อนช้าจะให้ตัวอย่างมีสหสัมพันธ์สูง ซึ่งให้ข้อมูลน้อยกว่าที่จำนวนระบุ หาก trace plot ดูเหมือนเลื่อนไปช้า ๆ ไม่ใช่แถบแนวนอนมีนอยส์ ให้ปรับจูนการกระจายข้อเสนอหรือใช้อัลกอริทึมที่ดีกว่า การเปลี่ยนจากเมโทรโปลิสไป NUTS มักแก้ได้ทันที
• ตัวอย่างน้อยเกินไป: ESS ต่ำหมายถึงค่าประมาณไม่น่าเชื่อถือ แม้โซ่จะบรรจบทางเทคนิคแล้ว หาก ESS เป็นเศษเสี้ยวของจำนวนดึงทั้งหมด ให้รันโซ่นานขึ้นหรือปรับปรุงการมิกซ์ก่อนสรุปผล

MCMC เทียบกับวิธีสุ่มตัวอย่างอื่น

MCMC ไม่ใช่วิธีเดียวในการประมาณการกระจาย ต่อไปนี้คือการเปรียบเทียบกับทางเลือกอื่น

• Rejection sampling ดึงข้อเสนอจากการกระจายที่ง่ายกว่าและยอมรับตามความสอดคล้องกับเป้าหมาย วิธีนี้ให้ผลแม่นยำเมื่อใช้ได้ แต่ในมิติสูงอัตราการยอมรับจะดิ่งลง
• Importance sampling ถ่วงน้ำหนักตัวอย่างจากการกระจายข้อเสนอเพื่อประมาณความคาดหมายภายใต้เป้าหมาย เร็วและเวิร์กในมิติต่ำ แต่ความแปรปรวนของน้ำหนักจะระเบิดเมื่อมิติสูง ในทางปฏิบัติ ตัวอย่างไม่กี่ตัวจะรับน้ำหนักเกือบทั้งหมด ทำให้ค่าประมาณไม่น่าเชื่อถือ
• Variational inference (VI) ฟิตการกระจายที่ง่ายกว่าเพื่อประมาณเป้าหมายโดยทำให้ความแตกต่างเชิงไดเวอร์เจนซ์ต่ำสุด VI เร็วกว่ามากและสเกลกับชุดข้อมูลใหญ่ได้ แต่มีข้อผิดพลาดจากการประมาณ การกระจายที่ฟิตอาจพลาดโครงสร้างของภายหลังที่ MCMC จะจับได้

สรุปสั้น ๆ คือ:

MCMC เทียบกับทางเลือกอื่น

MCMC เหมาะเมื่อความแม่นยำสำคัญกว่าความเร็ว หากต้องสเกลกับชุดข้อมูลใหญ่หรือทำอนุมานแบบเรียลไทม์ VI อาจคุ้มกับการแลกความแม่นยำ

สรุป

MCMC เป็นเครื่องมือที่ดูน่ากลัวจากภายนอก แต่เข้าใจง่ายขึ้นมากเมื่อเห็นว่ามันทำอะไรจริง ๆ นั่นคือสร้างสายโซ่ของตัวอย่างที่ค่อย ๆ สะท้อนรูปร่างของการกระจายที่คำนวณโดยตรงไม่ได้

และยิ่งเข้าใจง่ายเมื่อแยกเป็นส่วน: มาร์คอฟเชนกับวิธีมอนติคาร์โล

บทบาทของมันในสถิติเบย์เชียนยากจะพูดน้อยไปได้ การกระจายภายหลังที่เกินเอื้อมจะเข้าถึงได้ทันทีที่มีตัวสุ่มที่เชื่อถือได้ นั่นจึงเป็นเหตุผลที่ MCMC เป็นแกนกลางของไลบรารีโปรแกรมมิงเชิงความน่าจะเป็นอย่าง PyMC และ Stan

แต่ก่อนจะลงมือโค้ด ควรทำความเข้าใจสัญชาตญาณให้ถูกต้อง เข้าใจว่าทำไมโซ่ต้อง burn-in การมิกซ์คืออะไร และอ่าน trace plot อย่างไร โค้ดเองไม่ยาก เพราะไลบรารี Python ซ่อนความซับซ้อนทั้งหมดไว้หลังฟังก์ชันง่าย ๆ

หากต้องการชำนาญด้านการเรียนรู้ของเครื่อง ลงทะเบียนในเส้นทาง Machine Learning Scientist in Python ของเรา เนื้อหา 85 ชั่วโมงจะช่วยเตรียมความพร้อมสำหรับงานในปี 2026