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

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

โปรแกรมภาษาไพธอนนั้นเขียนง่าย อ่านง่าย แต่ขึ้นชื่อเรื่องความช้าเมื่อต้องทำงานมากๆเมื่อเทียบกับภาษาอื่นเช่น C หรือ Java

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

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

ถ้าโปรแกรมเราไม่ได้เขียนให้ใช้ Numpy แต่เราอยากให้มันวิ่งเร็วขึ้น โดยไม่อยากดัดแปลงโปรแกรม เราสามารถลองใช้ PyPy มาเรียกโปรแกรมของเรา เมื่อเราติดตั้ง PyPy แล้วเราสามารถเรียกโปรแกรมของเราชื่อ my_program.py ด้วยคำสั่ง pypy my_program.py หรือ pypy3 my_program.py แทน python my_program.py ได้เลย ส่วนใหญ่โปรแกรมจะเร็วขึ้นหลายเท่าเหมือนกัน

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

ยกตัวอย่างจากตอนที่ 1 และ 2 คือโปรแกรมทดลองแบ่งเส้นตรงเป็นสามส่วนแล้วดูว่าโอกาสที่ทั้งสามส่วนประกอบเป็นสามเหลี่ยมได้พอดีเป็นเท่าไร ถ้าเขียนด้วยไพธอนปกติก็มีหน้าตาแบบนี้:

import random

def prob_triangle(ntrials):
    """
    ประมาณความน่าจะเป็นที่เส้นตรงที่ถูกแบ่งเป็นสามส่วนแบบสุ่มๆ
    จะสามารถประกอบเป็นสามเหลี่ยมได้
    
    ทำการทดลอง ntrials ครั้ง
    """
    
    # ntriangles เก็บจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมได้
    ntriangles = 0 
    
    # ทำการทดลองแบ่งเส้นตรงความยาว 1 หน่วยเป็นสามส่วน
    # ทั้งหมดntrials ครั้ง
    
    for n in range(ntrials):
        
        # x, y คือตำแหน่งที่เราสุ่มตัดเส้นตรงสองตำแหน่ง
        # ทำให้แบ่งเส้นตรงเป็นสามส่วน
        # เราจะเรียงตำแหน่งให้ x <= y
        # ถ้าไม่เป็นอย่างนั้นเราจะสลับ x และy
        
        x = random.random()
        y = random.random()
        if x > y:
            x, y = y, x
        
        # เมื่อเราสุ่ม x, y มาตัดเส้นตรงได้แล้ว
        # เราจะมีเส้นตรงสามชิ้นยาว x, y-x, และ 1-y
        # เราจะเรียกชิ้นที่ยาวที่สุดว่า longest
        
        longest = max(x, y-x, 1-y)
        
        # ถ้าชิ้นที่ยาวที่สุดมีความยาวไม่เกิน 1/2 ของ
        # ความยาวเส้นตรงดั้งเดิม เราจะสามารถเอา
        # ทั้งสามชิ้นมาประกอบเป็นสามเหลี่ยมได้
        # และเราก็จะเพิ่มจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมสำเร็จ
        
        if longest < 0.5:
            ntriangles += 1
            
         # ประมาณความน่าจะเป็นที่เป็นสามเหลี่ยมสำเร็จ
         # เท่ากับ (จำนวนสามเหลี่ยม)/(จำนวนครั้งที่ทดลอง)
            
    return ntriangles/ntrials

ลองจับเวลาให้ทดลองกับเส้นตรง 10 ล้านเส้นใช้เวลา 3.67 วินาที:

%time prob_triangle(10_000_000)

CPU times: user 3.66 s, sys: 8.11 ms, total: 3.67 s
Wall time: 3.67 s

0.2499737

เราใช้ Numba โดยพิมพ์เพิ่ม 2 บรรทัดดังนี้ (บรรทัดที่ 1 และ 4)

from numba import njit #เรียกใช้ฟังก์ชั่นต่างๆใน numba
import random

@njit() #เพิ่มบรรทัดนี้เหนือฟังก์ชั่นที่เราเขียนไว้ มักจะทำงานได้เร็วขึ้น
def prob_triangle_numba(ntrials):
    """
    ประมาณความน่าจะเป็นที่เส้นตรงที่ถูกแบ่งเป็นสามส่วนแบบสุ่มๆ
    จะสามารถประกอบเป็นสามเหลี่ยมได้
    
    ทำการทดลอง ntrials ครั้ง
    """
    
    # ntriangles เก็บจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมได้
    ntriangles = 0 
    
    # ทำการทดลองแบ่งเส้นตรงความยาว 1 หน่วยเป็นสามส่วน
    # ทั้งหมดntrials ครั้ง
    
    for n in range(ntrials):
        
        # x, y คือตำแหน่งที่เราสุ่มตัดเส้นตรงสองตำแหน่ง
        # ทำให้แบ่งเส้นตรงเป็นสามส่วน
        # เราจะเรียงตำแหน่งให้ x <= y
        # ถ้าไม่เป็นอย่างนั้นเราจะสลับ x และy
        
        x = random.random()
        y = random.random()
        if x > y:
            x, y = y, x
        
        # เมื่อเราสุ่ม x, y มาตัดเส้นตรงได้แล้ว
        # เราจะมีเส้นตรงสามชิ้นยาว x, y-x, และ 1-y
        # เราจะเรียกชิ้นที่ยาวที่สุดว่า longest
        
        longest = max(x, y-x, 1-y)
        
        # ถ้าชิ้นที่ยาวที่สุดมีความยาวไม่เกิน 1/2 ของ
        # ความยาวเส้นตรงดั้งเดิม เราจะสามารถเอา
        # ทั้งสามชิ้นมาประกอบเป็นสามเหลี่ยมได้
        # และเราก็จะเพิ่มจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมสำเร็จ
        
        if longest < 0.5:
            ntriangles += 1
            
         # ประมาณความน่าจะเป็นที่เป็นสามเหลี่ยมสำเร็จ
         # เท่ากับ (จำนวนสามเหลี่ยม)/(จำนวนครั้งที่ทดลอง)
            
    return ntriangles/ntrials

เมื่อเราจับเวลาจะพบว่าฟังก์ชั่นแบบใช้ Numba ของเราเร็วขึ้นมาก ใช้เวลาเพียง 0.149 วินาที แทนที่จะเป็น 3.67 วินาที หรือเร็วเป็น 25 เท่า

%time prob_triangle_numba(10_000_000)

CPU times: user 148 ms, sys: 798 µs, total: 149 ms
Wall time: 149 ms

0.2499591

ถ้าเราเรียกฟังก์ชั่นเราครั้งเดียว ผลต่างระหว่าง 3 วินาที กับ 0.15 วินาทีจะดูไม่มาก แต่เมื่อเราให้คอมพิวเตอร์เรียกฟังก์ชั่นนั้นซ้ำๆกันเป็นพันครั้ง (เช่นเมื่อเราต้องการวาดฮิสโตแกรม) เวลาที่ใช้ก็จะต่างกันมากแบบ 1 ชั่วโมง vs. 2.5 นาที ดังนั้นเราควรรู้จักเทคนิคพวกนี้เมื่อสั่งให้คอมพิวเตอร์ทำงานหนักๆแทนเรา อาจจะประหยัดเวลาได้มาก

ถ้าคอมพิวเตอร์เรามีหลายๆคอร์ซีพียู (ซึ่งคอมปัจจุบันมักจะมีอย่างน้อยสองคอร์ขึ้นไปอยู่แล้ว) เราอาจขอให้ Numba พยายามใช้คอร์หลายๆคอร์ช่วยคำนวณด้วยคำสั่ง njit(parallel=True) และใช้ prange แทน range ในการคำนวณแบบนี้ (สังเกตบรรทัดที่ 1, 4, และ 23) เมื่อจับเวลาก็พบว่าเร็วขึ้นประมาณ 6 เท่าสำหรับซีพียูแบบ 8 คอร์:

from numba import njit, prange #เรียกใช้ฟังก์ชั่นต่างๆใน numba
import random

@njit(parallel=True) #เพิ่มบรรทัดนี้เหนือฟังก์ชั่นที่เราเขียนไว้ ฟังก์ชั่นเรามักจะทำงานได้เร็วขึ้น
def prob_triangle_numba_parallel(ntrials):
    """
    ประมาณความน่าจะเป็นที่เส้นตรงที่ถูกแบ่งเป็นสามส่วนแบบสุ่มๆ
    จะสามารถประกอบเป็นสามเหลี่ยมได้
    
    ทำการทดลอง ntrials ครั้ง
    """
    
    # ntriangles เก็บจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมได้
    ntriangles = 0 
    
    # ทำการทดลองแบ่งเส้นตรงความยาว 1 หน่วยเป็นสามส่วน
    # ทั้งหมดntrials ครั้ง
    
    # มีการใช้ prange แทน range เพื่อทำงานแบบพร้อมๆกัน
    # ตามจำนวนคอร์ใน CPU
    
    #ใช้ prange แทน range
    for n in prange(ntrials):
        
        # x, y คือตำแหน่งที่เราสุ่มตัดเส้นตรงสองตำแหน่ง
        # ทำให้แบ่งเส้นตรงเป็นสามส่วน
        # เราจะเรียงตำแหน่งให้ x <= y
        # ถ้าไม่เป็นอย่างนั้นเราจะสลับ x และy
        
        x = random.random()
        y = random.random()
        if x > y:
            x, y = y, x
        
        # เมื่อเราสุ่ม x, y มาตัดเส้นตรงได้แล้ว
        # เราจะมีเส้นตรงสามชิ้นยาว x, y-x, และ 1-y
        # เราจะเรียกชิ้นที่ยาวที่สุดว่า longest
        
        longest = max(x, y-x, 1-y)
        
        # ถ้าชิ้นที่ยาวที่สุดมีความยาวไม่เกิน 1/2 ของ
        # ความยาวเส้นตรงดั้งเดิม เราจะสามารถเอา
        # ทั้งสามชิ้นมาประกอบเป็นสามเหลี่ยมได้
        # และเราก็จะเพิ่มจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมสำเร็จ
        
        if longest < 0.5:
            ntriangles += 1
            
         # ประมาณความน่าจะเป็นที่เป็นสามเหลี่ยมสำเร็จ
         # เท่ากับ (จำนวนสามเหลี่ยม)/(จำนวนครั้งที่ทดลอง)
            
    return ntriangles/ntrials

ถ้าจะใช้ parallel และ prange ควรอ่านเอกสารคู่มือเรื่องนี้ของ Numba ก่อนนะครับว่ามีข้อควรระวังอะไรบ้าง

หวังว่าผู้อ่านจะได้ไอเดียหรือประโยชน์บ้าง ถ้ามีคำแนะนำหรือข้อสงสัยส่งข้อความอินบ๊อกซ์ไปที่เพจวิทย์พ่อโก้บนเฟซบุ๊คได้ครับ: https://www.facebook.com/witpokosci/

(ตอนที่แล้วอยู่ที่ https://witpoko.com/?p=7485)

เราคุยกันต่อจากตอนที่แล้วที่เราให้คอมพิวเตอร์ทำการทดลองหาคำตอบที่เราไม่รู้ (computer simulation)

คำถามคือเอาเส้นตรงมาหนึ่งเส้น แบ่งเป็นสามท่อน ความน่าจะเป็นที่ทั้งสามท่อนสามารถมาประกอบกันเป็นสามเหลี่ยมได้พอดีมีเท่าไร คำตอบจริงๆจากทฤษฎีคือ 1/4 แต่เนื่องจากเราไม่รู้คำตอบจากทฤษฎีเราจึงให้คอมพิวเตอร์ทดลองแบ่งเส้นตรงเป็นสามท่อนแบบสุ่มๆหลายๆครั้งแล้วนับจำนวนครั้งที่ประกอบเป็นสามเหลี่ยมได้สำเร็จ เปรียบเทียบกับจำนวนการทดลองทั้งหมด อัตราส่วนก็คือค่าความน่าจะเป็นโดยประมาณ

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

(ถึงตอนนี้ถ้าอยู่ในห้องเรียน เราก็ควรถามเด็กๆว่าเด็กๆคิดว่าควรทำยังไงดี ให้เขาคิดและเสนอความคิดกัน)

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

nexperiments = 10
ntrials = 100_000
print(f"สำหรับการทดลองแต่ละการทดลอง")
print(f"จะแบ่งเส้นตรงแบบสุ่มๆเป็นสามชิ้นทั้งหมด {ntrials:,} ครั้ง")
print("แล้วดูว่าทั้งสามชิ้นสามารถประกอบเป็นสามเหลี่ยมได้ไหม\n")
for e in range(1,nexperiments+1):
    print(f"การทดลองครั้งที่ {e}, ความน่าจะเป็นที่ได้สามเหลี่ยม = {prob_triangle(ntrials):.3f}")

ผลที่ได้คือ:

สำหรับการทดลองแต่ละการทดลอง
จะแบ่งเส้นตรงแบบสุ่มๆเป็นสามชิ้นทั้งหมด 100,000 ครั้ง
แล้วดูว่าทั้งสามชิ้นสามารถประกอบเป็นสามเหลี่ยมได้ไหม

การทดลองครั้งที่ 1, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.249
การทดลองครั้งที่ 2, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.249
การทดลองครั้งที่ 3, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.250
การทดลองครั้งที่ 4, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.250
การทดลองครั้งที่ 5, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.253
การทดลองครั้งที่ 6, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.252
การทดลองครั้งที่ 7, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.250
การทดลองครั้งที่ 8, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.251
การทดลองครั้งที่ 9, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.246
การทดลองครั้งที่ 10, ความน่าจะเป็นที่ได้สามเหลี่ยม = 0.252

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

วิธีดูว่าหน้าตาคำตอบเป็นอย่างไรที่มีประโยชน์มากคือดูด้วยฮิสโตแกรม (histogram) ซึ่งจะบอกว่าค่าคำตอบที่ได้แต่ละค่าเกิดบ่อยแค่ไหน

ก่อนอื่นเราต้องทำการทดลองหลายๆครั้งก่อน เช่นทำสิบครั้ง:

nexperiments = 10 #จำนวนครั้งที่ทำการทดลอง
ntrials = 100_000 #จำนวนเส้นตรงที่ลองแบ่งเป็นสามส่วนในการทดลองแต่ละครั้ง
answers = [] #เก็บคำตอบจากการทดลองทั้งหมดไว้ในนี้
for e in range(1,nexperiments+1):
    answers.append(prob_triangle(ntrials))

print(answers)

คำตอบเราเก็บไว้ในลิสต์ชื่อ answers ซึ่งเราสามารถเอาไปทำอะไรต่อได้เช่นวาดฮิสโตแกรม ตอนนี้เราทำการทดลองแค่ 10 ครั้งมีหน้าตาประมาณนี้:

[0.2516, 0.251, 0.25219, 0.25302, 0.25189, 0.25049, 0.24938, 0.24922, 0.24566, 0.2502]

เราเอามาวาดฮิสโตแกรมได้ด้วย matplotlib แบบนี้ (หรือใช้อย่างอื่นเช่น Excel วาดก็ได้):

%matplotlib inline
import matplotlib.pyplot as plt #เรียกใช้ฟังก์ชั่นวาดกราฟต่างๆ

plt.hist(answers) #เอาข้อมูล answers ใส่ฮิสโตแกรม
plt.xlabel("ความน่าจะเป็นที่เป็นสามเหลี่ยม", fontname="Tahoma")
plt.ylabel("จำนวนครั้งที่เกิด",fontname="Tahoma")
plt.show() #วาดฮิสโตแกรม

ฮิสโตแกรมที่ได้หน้าตาตลกดี เพราะจำนวนครั้งที่ทดลองมีเพียงสิบครั้ง แต่อย่างไรก็ตามคำตอบทั้งสิบอยู่ในช่วง 0.246 ถึง 0.253:

เราเพิ่มจำนวนครั้งที่ทดลองเป็น 100 ครั้งโดยเปลี่ยนตัวแปร nexperiments จาก 10 เป็น 100 แล้ววาดฮิสโตแกรมใหม่:

nexperiments = 100 #จำนวนครั้งที่ทำการทดลอง
ntrials = 100_000 #จำนวนเส้นตรงที่ลองแบ่งเป็นสามส่วนในการทดลองแต่ละครั้ง
answers = [] #เก็บคำตอบจากการทดลองทั้งหมดไว้ในนี้
for e in range(1,nexperiments+1):
    answers.append(prob_triangle(ntrials))

plt.hist(answers) #เอาข้อมูล answers ใส่ฮิสโตแกรม
plt.xlabel("ความน่าจะเป็นที่เป็นสามเหลี่ยม", fontname="Tahoma")
plt.ylabel("จำนวนครั้งที่เกิด",fontname="Tahoma")
plt.show() #วาดฮิสโตแกรม

หน้าตาฮิสโตแกรมดูดีขี้น ค่าตรงกลางๆสูงกว่าขอบๆอย่างเห็นได้ชัด

ถ้าเราทำการทดลองหลายๆครั้งมากๆ ให้เป็น 10,000 ครั้งแทนที่จะเป็น 100 ครั้ง ฮิสโตแกรมเราจะหน้าตาแบบนี้ ค่าคำตอบกระจุกตัวอยู่แถวๆ 0.25 ไม่มีตำตอบไหนอยู่ห่างไกลนัก:

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

หลังจากเราดูภาพฮิสโตแกรมว่าคำตอบอยู่เป็นกลุ่ม ไม่กระจัดกระจายเกินไป เราอาจจะสรุปคำตอบที่ได้มาด้วยค่าเฉลี่ย (mean) และค่าเบี่ยงเบนมาตรฐาน (standard deviation) เช่นการคำนวณด้วย Numpy แบบนี้ครับ:

import numpy as np #เรียกใช้ฟังก์ชั่นต่างๆใน numpy
print(f"จากการทดลอง {nexperiments:,} ครั้ง")
print(f"ค่าเฉลี่ยของคำตอบคือ {np.mean(answers):.3f}")
print(f"ค่าเบี่ยงเบนมาตรฐานคือ {np.std(answers):.3f}")

จากการทดลอง 1,000 ครั้ง
ค่าเฉลี่ยของคำตอบคือ 0.250
ค่าเบี่ยงเบนมาตรฐานคือ 0.001

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

หวังว่าผู้อ่านจะได้ไอเดียหรือประโยชน์บ้าง ถ้ามีคำแนะนำหรือข้อสงสัยส่งข้อความอินบ๊อกซ์ไปที่เพจวิทย์พ่อโก้บนเฟซบุ๊คได้ครับ: https://www.facebook.com/witpokosci/

(ตอนที่แล้วอยู่ที่ https://witpoko.com/?p=7485 ตอนต่อไปอยู่ที่ https://witpoko.com/?p=7526)