গাণিতিক কাজ, aggregation আর random — ufunc, axis আর seed

আগের তিনটি পর্বে (পর্ব ১৩-এ array, shape আর dtype; পর্ব ১৪-এ indexing আর slicing; পর্ব ১৫-এ vectorization আর broadcasting) দেখেছেন NumPy array কী আর কীভাবে তার ভেতরের অংশ বেছে নেওয়া যায়। এই পর্বে সেই array দিয়ে আসল হিসাব: গাণিতিক function, সারাংশ (মোট, গড়, ছড়ানো) আর random সংখ্যা। শেষে এই তিনটি মিলিয়ে বাজারঘরের (এই কোর্সের কাল্পনিক অনলাইন দোকান) একটি ব্যবসায়িক প্রশ্নের উত্তর খুঁজবেন: একটি অর্ডার ৩ দিনের বেশি সময় নেওয়ার সম্ভাবনা কত? সব কোড পাতাতেই চলে; প্রথমবার চালালে NumPy নামতে কয়েক সেকেন্ড লাগতে পারে। Seed-এর মৌলিক ধারণা আগেই এসেছে পর্ব ১১-এ; এখানে তার পরের ধাপ।

১ · Ufunc — পুরো array-এ একসঙ্গে গাণিতিক কাজ

বাজারঘরের গুদামে ৫০০টি পণ্যের দাম একটি array-তে আছে। সব দামে ১৫% ভ্যাট যোগ করতে হবে, ঋণাত্মক ভুল দাম শূন্য করতে হবে, আর বিশ্লেষণের জন্য দামের logarithm লাগবে। খাঁটি Python-এ প্রতিটি কাজের জন্য একটি loop লাগত। NumPy-তে এক লাইন: np.log(price). এই ধরনের function-এর নাম universal function, ছোট করে ufunc: এমন একটি function, যা একটি array নেয় আর প্রতিটি উপাদানের ওপর আলাদাভাবে একই কাজ করে একই shape-এর নতুন array ফেরত দেয়। Loop-টি চলে NumPy-র ভেতরের compiled C কোডে, যা পর্ব ৯-এ মেপে দেখা vectorization-এর মূল কারণ।

সবচেয়ে বেশি কাজে লাগে এমন ufunc-গুলো এক নজরে:

খেয়াল করুন, np.maximum আর np.max এক জিনিস নয়। প্রথমটি দুটি array-কে ঘরে ঘরে তুলনা করে একই আকারের array দেয়; দ্বিতীয়টি একটি array-এর সব সংখ্যার মধ্যে সবচেয়ে বড়টি, মানে একটি সারাংশ (অংশ ৩)। +, *, >-এর মতো চিহ্নও আসলে ufunc: a + b লিখলে ভেতরে np.add(a, b) চলে। দুটি array-এর shape আলাদা হলে broadcasting-এর নিয়ম (পর্ব ১৫) ঠিক করে কোন ঘর কার সঙ্গে মিলবে।

out= আর where=

প্রতিটি ufunc দুটি কাজের argument নেয়। out= বলে দেয় ফল কোন আগে-বানানো array-তে লেখা হবে; বড় array-তে বারবার নতুন মেমরি না নিয়ে একই জায়গা ব্যবহারে এটি কাজে লাগে। where= একটি boolean mask (পর্ব ১৪) নেয়: যেখানে True, শুধু সেখানে কাজ হয়; যেখানে False, সেখানে out-এ আগের মানই থাকে। তাই দুটি একসঙ্গে দিলে শূন্য দিয়ে ভাগ এড়ানোর একটি পরিচ্ছন্ন উপায় পাওয়া যায়: প্রতিটি পণ্যের “বিক্রির হার” = বিক্রি ÷ মজুত, কিন্তু মজুত শূন্য হলে হার ০ ধরা হবে।

নিচের কোডটি লেখকের যন্ত্রে সত্যিকারের run (Python 3.14.4, NumPy 2.5.1); এখানে ব্রাউজারের ব্লক না রেখে কাগজে দেওয়া হলো, কারণ পাতার চারটি চালানো যায় এমন ব্লক পরের অংশগুলোর জন্য রাখা। আপনি কোডটি পরের যেকোনো ব্লকে বসিয়ে চালিয়ে মেলাতে পারেন।

~/bazarghar-ai$ python ufunc_demo.py

import numpy as np
price = np.array([120.0, 85.5, 300.0, 49.99])        # টাকা
print(np.sqrt(np.array([4.0, 9.0, 2.0])))
print(np.exp(np.array([0.0, 1.0, 2.0])))
print(np.log(np.array([1.0, np.e, 100.0])))
print(np.abs(np.array([-3, 0, 5])))
print(np.round(np.array([2.5, 3.5, 0.125]), 2), np.round(np.array([2.5, 3.5])))
print(np.clip(np.array([-2, 3, 12, 7]), 0, 10))
print(np.maximum(np.array([3, 8, 1]), np.array([5, 2, 4])))
print(type(np.sqrt), np.add.nin, np.add.nout)

out = np.empty(4)
r = np.multiply(price, 1.15, out=out)                 # ভ্যাটসহ দাম
print(r is out, out)

stock = np.array([10.0, 0.0, 4.0, 0.0])
sold = np.array([5.0, 0.0, 1.0, 3.0])
rate = np.zeros(4)
np.divide(sold, stock, out=rate, where=stock > 0)
print(rate)

[2.         3.         1.41421356]
[1.         2.71828183 7.3890561 ]
[0.         1.         4.60517019]
[3 0 5]
[2.5  3.5  0.12] [2. 4.]
[ 0  3 10  7]
[5 8 4]
<class 'numpy.ufunc'> 2 1
True [138.      98.325  345.      57.4885]
[0.5  0.   0.25 0.  ]

প্রতিটি দুই-input ufunc-এর সঙ্গে কয়েকটি method-ও আসে, যা পরের অংশের সেতু: np.add.reduce(a) সব উপাদান যোগ করে একটি সংখ্যা দেয় (এটিই np.sum), আর np.add.accumulate(a) চলতি যোগফল দেয় (এটিই np.cumsum). মানে সারাংশের function-গুলো আসলে ufunc-কে বারবার প্রয়োগ করা।

২ · Float-এর কিনারা — inf, nan আর শূন্য দিয়ে ভাগ

Python-এ 1 / 0 লিখলে প্রোগ্রাম থেমে যায় ZeroDivisionError দিয়ে। NumPy-র float array-তে একই ভাগ থামে না; ফলের ঘরে একটি বিশেষ মান বসে, আর সাধারণত একটি সতর্কবার্তা (warning) ছাপা হয়। কারণ NumPy floating-point সংখ্যার আন্তর্জাতিক নিয়ম IEEE 754 মেনে চলে, যেখানে এমন কিছু “সংখ্যা” আছে যা সাধারণ সংখ্যা নয়:

• inf (অসীম): 900.0 / 0.0 হলো inf, -900.0 / 0.0 হলো -inf. খুব বড় ফল, যেমন np.exp(1000.0), float-এ না ধরলে সেটিও inf হয় (একে বলে overflow).
• nan (Not a Number, সংখ্যা নয়): যার কোনো অর্থপূর্ণ উত্তর নেই, যেমন 0.0 / 0.0, inf - inf, বা np.sqrt(-1.0). ডেটায় “ফাঁকা ঘর” বোঝাতেও nan ব্যবহার হয়; pandas-এ CSV-র খালি ঘর পড়লে প্রায়ই nan আসে (পর্ব ২৩)।

nan-এর সবচেয়ে বড় ফাঁদ: nan কোনো কিছুর সমান নয়, নিজেরও নয়। np.nan == np.nan হলো False. তাই if x == np.nan: কখনো সত্যি হয় না, আর ফাঁকা ঘর খুঁজতে এভাবে লিখলে কোড চুপচাপ কিছুই খুঁজে পায় না। সঠিক পরীক্ষা np.isnan(x); আর np.isfinite(x) একসঙ্গে nan আর দুই দিকের inf বাদ দেয়। দ্বিতীয় ফাঁদ: nan ছোঁয়াচে। যেকোনো হিসাবে একটি nan ঢুকলে ফলও nan: nan + 5, nan * 0, আর পুরো array-এর sum-ও।

সতর্কবার্তা আর np.errstate

শূন্য দিয়ে ভাগে NumPy যে সতর্কবার্তা দেয়, তা লেখকের যন্ত্রে সাধারণ Python-এ চালালে দেখতে এমন (সত্যিকারের run, NumPy 2.5.1):

<string>:3: RuntimeWarning: divide by zero encountered in divide
<string>:3: RuntimeWarning: invalid value encountered in divide
<string>:4: RuntimeWarning: divide by zero encountered in log
<string>:4: RuntimeWarning: invalid value encountered in sqrt

“divide by zero” মানে শূন্য দিয়ে ভাগে inf এসেছে; “invalid value” মানে nan এসেছে। NumPy-র উত্তর np.errstate, একটি context manager (Python-এর with ব্লক; বিস্তারিত Python কোর্সের context manager পর্বে)। এর ভেতরে প্রতিটি ধরনের সমস্যার জন্য আচরণ ঠিক করা যায়: "ignore" (চুপ থাকা), "warn" (সতর্ক করা, default), "raise" (FloatingPointError দিয়ে থেমে যাওয়া)। ধরনগুলো হলো divide, invalid, over আর under.

Float তুলনা — np.isclose

পর্ব ৯-এ দেখেছেন, 0.1 + 0.2 == 0.3 হলো False, আর Python-এ তুলনা হয় math.isclose দিয়ে। NumPy-র সমতুল্য np.isclose(a, b), যা ঘরে ঘরে কাজ করে, আর পুরো array একসঙ্গে মেলাতে np.allclose. শর্তটি হলো |a − b| ≤ atol + rtol × |b|, যেখানে rtol হলো আপেক্ষিক সহনশীলতা (default 1e-05, মানে b-এর এক লাখ ভাগের এক ভাগ) আর atol হলো পরম সহনশীলতা (default 1e-08). একটি সূক্ষ্ম ফাঁদ: atol-এর কারণে খুব ছোট সংখ্যা শূন্যের “কাছাকাছি” ধরা পড়ে। np.isclose(1e-9, 0.0) হলো True; আপনার ডেটা যদি নিজেই খুব ছোট এককে হয় (যেমন সম্ভাবনা 10^-9), তাহলে atol=0 দিন। চালানো যায় এমন উদাহরণ অংশ ৫-এ।

৩ · Aggregation — অনেক সংখ্যা থেকে একটি সারাংশ

ক্রিকেটের স্কোরকার্ডের কথা ভাবুন। প্রতিটি বলের রান লেখা থাকে, কিন্তু খেলা শেষে আমরা জানতে চাই মোট রান, গড়, সর্বোচ্চ ইনিংস। অনেক সংখ্যা থেকে একটি সংখ্যায় আসার এই কাজের নাম aggregation. ufunc array-এর shape বদলায় না; aggregation এক বা একাধিক দিক ভেঙে ছোট করে। একটি উদাহরণ দিয়ে সবগুলো হাতে হিসাব করি: ঢাকা থেকে চার দিনে (শনি থেকে মঙ্গল) বাজারঘরের অর্ডার এসেছে x = [12, 15, 9, 14].

প্রতিটি সংখ্যা লেখকের যন্ত্রে NumPy দিয়ে মেলানো। এবার শব্দগুলোর মানে। Median হলো সাজানো তালিকার ঠিক মাঝের মান; জোড় সংখ্যক মান থাকলে মাঝের দুটির গড়। Variance (ভেদাঙ্ক, চিহ্ন σ²) হলো গড় থেকে দূরত্বের বর্গের গড়, আর standard deviation (আদর্শ বিচ্যুতি, σ) তার বর্গমূল, যা মূল এককে (এখানে “অর্ডার”) ছড়ানো মাপে। Percentile p মানে এমন একটি মান, যার নিচে মোটামুটি p% ডেটা; ২৫, ৫০ আর ৭৫ percentile-কে quartile বলে, আর ৫০ percentile-ই median. NumPy-র default পদ্ধতি দুই পাশের মানের মাঝে সরলরেখায় ভাগ করে (method="linear"), যা ওপরের হিসাবে দেখলেন; অন্য সফটওয়্যার অন্য পদ্ধতি নিতে পারে, তাই ছোট ডেটায় সংখ্যা সামান্য আলাদা আসতে পারে। cumsum (cumulative sum) চলতি যোগফল: মাসের প্রথম দিন থেকে আজ পর্যন্ত মোট বিক্রি।

ddof — ৪ দিয়ে ভাগ, নাকি ৩ দিয়ে?

ddof (delta degrees of freedom) বলে দেয় ভাগ হবে n − ddof দিয়ে। NumPy-র default ddof=0: চারটি সংখ্যাই যদি আপনার পুরো জগৎ (population) হয়, তাহলে n দিয়ে ভাগ ঠিক। কিন্তু এই চার দিন যদি সারা বছরের একটি নমুনা (sample) হয়, আর আপনি সারা বছরের ছড়ানো অনুমান করতে চান, তাহলে n − 1 দিয়ে ভাগ করলে অনুমানটি গড়ে পক্ষপাতহীন হয় (variance-এর জন্য; একে Bessel's correction বলে)। ব্যবহারিক ফাঁদটি হলো: pandas-এর std আর var-এর default ddof=1, NumPy-র 0. একই ডেটায় দুই library দুই উত্তর দেয় (2.291 বনাম 2.646), আর কেউ ভুল নয়। নমুনা, পক্ষপাতহীন অনুমান আর এর তত্ত্ব আসবে পর্ব ৩৭-এ।

৪ · axis — জেলা × দিনের টেবিলে কোন দিকে হিসাব

এক জেলার চার দিনের ডেটা ছিল এক মাত্রার। আসল ডেটা প্রায় সবসময় টেবিল: বাজারঘরের তিনটি জেলা (ঢাকা, চট্টগ্রাম, খুলনা) আর চারটি দিন (শনি, রবি, সোম, মঙ্গল)। সারি জেলা, কলাম দিন, তাই shape (3, 4). এখন orders.sum() লিখলে সব ১২টি ঘর যোগ হয়ে একটি সংখ্যা আসে। কিন্তু ম্যানেজার জানতে চান দুটি আলাদা প্রশ্ন: “প্রতিটি দিনে মোট কত?” আর “প্রতিটি জেলায় মোট কত?” এটিই axis argument-এর কাজ।

নিয়মটি একটি বাক্যে: axis=k মানে k নম্বর মাত্রাটি ভেঙে যাবে (collapse), বাকিগুলো থাকবে। shape (3, 4)-এ মাত্রা ০ হলো সারি (৩টি জেলা), মাত্রা ১ হলো কলাম (৪টি দিন)।

• axis=0: সারির দিকটি ভাঙে। হিসাব চলে ওপর থেকে নিচে, প্রতিটি কলামের ভেতরে। ফল প্রতিটি দিনের জন্য একটি সংখ্যা, shape (4,).
• axis=1: কলামের দিকটি ভাঙে। হিসাব চলে বাঁ থেকে ডানে, প্রতিটি সারির ভেতরে। ফল প্রতিটি জেলার জন্য একটি সংখ্যা, shape (3,).

অনেকে উল্টো মনে রাখেন: “axis=0 মানে সারি, তাই সারির যোগফল।” না; axis=0 মানে সারিগুলো মিশে যায়, তাই যা বাকি থাকে তা কলাম। মনে রাখার সহজ উপায়: ফলের shape লিখুন, আর axis-এর ঘরটি কেটে দিন। (3, 4) থেকে axis=0 কাটলে (4,). নিচের ধাপে দেখুন।

এবার একই হিসাব সত্যিই চালিয়ে দেখুন, প্রতিটি ফলের shape সহ। সঙ্গে mean(axis=1) (প্রতিটি জেলার দৈনিক গড়) আর argmax (সবচেয়ে বড়টির index) দিয়ে “সেরা দিন” খোঁজা। শেষের দুই লাইনে keepdims, যা পরের অংশে ধাপে ধাপে দেখবেন।

import numpy as np

districts = ["ঢাকা", "চট্টগ্রাম", "খুলনা"]
days = ["শনি", "রবি", "সোম", "মঙ্গল"]
orders = np.array([[12, 15,  9, 14],      # সারি = জেলা
                   [ 8,  6, 10,  7],      # কলাম = দিন
                   [ 5,  7,  6,  3]], dtype=np.int64)
print("table shape:", orders.shape)

by_day = orders.sum(axis=0)               # ওপর থেকে নিচে, প্রতিটি কলাম
by_district = orders.sum(axis=1)          # বাঁ থেকে ডানে, প্রতিটি সারি
print("axis=0:", by_day, by_day.shape)
print("axis=1:", by_district, by_district.shape)
print("সব মিলিয়ে:", orders.sum(), "| দুই পথে:", by_day.sum(), by_district.sum())

for name, avg in zip(districts, orders.mean(axis=1)):
    print(f"{name} দিনে গড় {avg:.2f}টি অর্ডার")
print("সেরা দিন:", days[orders.sum(axis=0).argmax()])
print("প্রতিটি জেলার সেরা দিন:", [days[i] for i in orders.argmax(axis=1)])

tot = orders.sum(axis=1, keepdims=True)   # shape (3, 1)
print("keepdims:", tot.shape, "→ ভাগের পর", (orders / tot).shape)
np.set_printoptions(precision=3)
print(orders / tot)
print("প্রতি সারির যোগ:", (orders / tot).sum(axis=1))

table shape: (3, 4)
axis=0: [25 28 25 24] (4,)
axis=1: [50 31 21] (3,)
সব মিলিয়ে: 102 | দুই পথে: 102 102
ঢাকা দিনে গড় 12.50টি অর্ডার
চট্টগ্রাম দিনে গড় 7.75টি অর্ডার
খুলনা দিনে গড় 5.25টি অর্ডার
সেরা দিন: রবি
প্রতিটি জেলার সেরা দিন: ['রবি', 'সোম', 'রবি']
keepdims: (3, 1) → ভাগের পর (3, 4)
[[0.24  0.3   0.18  0.28 ]
 [0.258 0.194 0.323 0.226]
 [0.238 0.333 0.286 0.143]]
প্রতি সারির যোগ: [1. 1. 1.]

তিন বা তার বেশি মাত্রাতেও নিয়ম একই: (7, 3, 4) (সপ্তাহ × জেলা × দিন)-এ sum(axis=(0, 2)) প্রতিটি জেলার মোট দেয়, shape (3,), আর axis=-1 মানে শেষ মাত্রা। পর্ব ২৪-এর groupby-এও একই ভাবনা: কোন দিক মিশবে, কোন দিক থাকবে।

৫ · keepdims আর ফাঁকা দিন — shape ধরে রাখা আর nan-সচেতন হিসাব

keepdims=True কেন লাগে

ওপরের কোডের শেষ অংশে প্রশ্নটি ছিল: প্রতিটি জেলার নিজের মোটের কত ভাগ কোন দিনে? মানে টেবিলের প্রতিটি সারিকে সেই সারির যোগফল দিয়ে ভাগ। স্বাভাবিকভাবে লিখতে ইচ্ছা করে orders / orders.sum(axis=1). কিন্তু এটি error দেয়, কারণ broadcasting-এর নিয়ম shape-গুলো ডান দিক থেকে মেলায় (পর্ব ১৫)। keepdims=True দিলে ভাঙা মাত্রাটি মুছে না গিয়ে দৈর্ঘ্য ১ হয়ে থেকে যায়, আর তখন broadcasting ঠিক যেমন চাই তেমন কাজ করে। ধাপে ধাপে দেখুন।

1tot = orders.sum(axis=1)2orders / tot          # (3,4) / (3,)3tot = orders.sum(axis=1, keepdims=True)4share = orders / tot  # (3,4) / (3,1)5share.sum(axis=1)

একটি বিপজ্জনক ব্যতিক্রম মনে রাখুন: টেবিলটি যদি বর্গাকার হতো, ধরুন (4, 4), তাহলে orders / orders.sum(axis=1) কোনো error দিত না। (4, 4) আর (4,) ডান থেকে মেলে, তাই NumPy চুপচাপ প্রতিটি কলামকে ভাগ করত একটি সারির মোট দিয়ে। ফল ভুল, কিন্তু চোখে পড়ার মতো কোনো লক্ষণ নেই। এ কারণেই কোনো সারাংশ দিয়ে আবার মূল টেবিলে হিসাব করতে হলে keepdims=True অভ্যাস করুন, আর শেষে একটি যাচাই (যেমন সারির যোগ 1) রাখুন।

ফাঁকা দিন — nanmean আর বন্ধুরা

বাজারঘরের একটি শাখার পাঁচ দিনের বিক্রির খাতা থেকে মঙ্গলবারের পাতাটি হারিয়ে গেছে। ডেটায় সেই ঘরে nan. এখন sales.mean() দেয় nan, কারণ nan ছোঁয়াচে। NumPy-র সমাধান nan-সচেতন function-এর একটি পরিবার: np.nanmean, np.nansum, np.nanmedian, np.nanstd, np.nanmin, np.nanmax, np.nanpercentile. এরা nan-গুলোকে বাদ দিয়ে বাকিগুলোর হিসাব করে। কিন্তু “বাদ দেওয়া” একটি সিদ্ধান্ত: হারানো দিনটি হয়তো ঈদের ছুটি ছিল (বিক্রি আসলে শূন্য), বা হয়তো রেকর্ড হয়নি (বিক্রি হয়েছে কিন্তু জানা নেই)। প্রথম ক্ষেত্রে nan নয়, 0 বসানো উচিত। ফাঁকা ঘর পূরণ করার কৌশল (imputation) pandas দিয়ে পর্ব ২৩-এ।

নিচের কোডে চারটি জিনিস: হারানো দিন আর nanmean; nan == nan-এর ফাঁদ; শূন্য দিয়ে ভাগ, প্রথমে np.errstate দিয়ে ইচ্ছা করে চুপ করানো, তারপর where= দিয়ে এড়ানো, তারপর "raise" দিয়ে ধরার চেষ্টা (এর একটি চমক আছে); আর np.isclose.

import numpy as np

sales = np.array([4200.0, 3900.0, np.nan, 5100.0, 4800.0])  # মঙ্গলবার খাতা হারিয়েছে
print("mean:", sales.mean(), "| nanmean:", np.nanmean(sales))
print("কয়টি ফাঁকা:", np.isnan(sales).sum())
print("nan == nan ?", np.nan == np.nan, "| sales[2] == np.nan ?", sales[2] == np.nan)
print("ঠিক পরীক্ষা:", np.isnan(sales[2]))

revenue = np.array([0.0, 900.0, 1500.0])
orders = np.array([0.0, 0.0, 6.0])
with np.errstate(divide="ignore", invalid="ignore"):
    per_order = revenue / orders            # 0/0 আর 900/0
print("প্রতি অর্ডার:", per_order)
print("finite?", np.isfinite(per_order))
safe = np.divide(revenue, orders, out=np.zeros(3), where=orders > 0)
print("where দিয়ে নিরাপদ:", safe)

try:
    with np.errstate(divide="raise"):
        revenue[1] / orders[1:2]
except FloatingPointError as e:
    print("ধরা পড়ল:", e)                    # সাধারণ Python-এ
else:
    print("কিছু ধরা পড়েনি (ব্রাউজারের Python)")

a = 0.1 + 0.2
print(a == 0.3, np.isclose(a, 0.3), np.isclose(1e-9, 0.0), np.isclose(1e-9, 0.0, atol=0))
print("inf + 1 =", np.inf + 1, "| inf - inf =", np.inf - np.inf)

mean: nan | nanmean: 4500.0
কয়টি ফাঁকা: 1
nan == nan ? False | sales[2] == np.nan ? False
ঠিক পরীক্ষা: True
প্রতি অর্ডার: [ nan  inf 250.]
finite? [False False  True]
where দিয়ে নিরাপদ: [  0.   0. 250.]
কিছু ধরা পড়েনি (ব্রাউজারের Python)
False True True False
inf + 1 = inf | inf - inf = nan

৬ · Random-এর নতুন পথ — default_rng আর Generator-এর কাজ

পর্ব ১১-এ দেখেছেন, কম্পিউটারের random সংখ্যা আসে একটি PRNG থেকে, যা একটি seed থেকে শুরু করে একটি নির্দিষ্ট ধারা দেয়। এখানে দেখব NumPy-তে সেই যন্ত্রটি কীভাবে সাজানো, আর ডেটার কাজে কোন method কখন লাগে।

NumPy-র আধুনিক random ব্যবস্থার তিনটি স্তর। নিচের ছবিতে ওপর থেকে নিচে পড়ুন।

কেন পুরোনো np.random.seed নয়

অনেক পুরোনো টিউটোরিয়াল আর AI সহকারীর কোডে এখনো দেখবেন np.random.seed(42) আর np.random.randint(...). এগুলো legacy (পুরোনো, শুধু পেছনের সামঞ্জস্যের জন্য রাখা) পদ্ধতি: পুরো প্রোগ্রামে একটিমাত্র লুকানো global অবস্থা, যার ভেতরে Mersenne Twister (MT19937) নামের PRNG. সমস্যা তিনটি। এক, কোনো library-র ভেতরের একটি ডাক সেই একই ধারা থেকে সংখ্যা নিলে আপনার পরের সব সংখ্যা সরে যায়, আর কেন, তা কোড দেখে বোঝা যায় না। দুই, দুটি আলাদা কাজের জন্য আলাদা ধারা রাখা কঠিন। ২০১৯ সালে NumPy 1.17-এ নতুন ব্যবস্থা আসে (NumPy-র নকশা দলিল NEP 19-এ এর যুক্তি আছে): default_rng একটি নিজস্ব Generator ফেরত দেয়, যা আপনি একটি সাধারণ variable-এর মতো function-এ পাঠাতে পারেন। Default BitGenerator হলো PCG64, যা ছোট অবস্থায় দ্রুত আর পরিসংখ্যানের পরীক্ষায় ভালো ফল দেয় বলে বেছে নেওয়া হয়েছিল। NumPy-র docs নতুন কোডে Generator ব্যবহারের পরামর্শ দেয়।

ডেটার কাজের প্রধান method

choice-এর দুটি রূপের পার্থক্য বাস্তবে গুরুত্বপূর্ণ। ফেরত দিয়ে (with replacement, default) মানে লটারির বাক্স থেকে একটি কাগজ তুলে নাম লিখে আবার বাক্সে রাখা: একই জিনিস বারবার আসতে পারে। ফেরত ছাড়া (replace=False) মানে কাগজ আর ফেরত যায় না: ৫ জন কর্মীর মধ্যে ৩ জনকে আলাদা আলাদা কাজে বাছাই, বা ডেটা থেকে train আর test সারি বাছাই, যেখানে একই সারি দুবার আসা চলবে না। আর p= দিয়ে weighted sampling: বাজারঘরে যদি ৬০% অর্ডার ক্যাশ অন ডেলিভারি (COD), ৩০% bKash আর ১০% কার্ডে হয়, তাহলে p=[0.6, 0.3, 0.1] দিলে বানানো অর্ডারগুলোতেও প্রায় সেই অনুপাত আসে।

import numpy as np

rng_a = np.random.default_rng(42)
rng_b = np.random.default_rng(42)          # একই seed, আলাদা generator
print("একই seed:", rng_a.integers(1, 7, size=5), rng_b.integers(1, 7, size=5))
print("পরের ৫টি:", rng_a.integers(1, 7, size=5))   # stream এগিয়ে গেছে

rng = np.random.default_rng(2026)
print("random():", rng.random(3).round(3))         # [0, 1)
demand = rng.normal(loc=40, scale=6, size=7)       # দৈনিক চাহিদা
print("normal:", demand.round(1))
print("দৈনিক অর্ডার (পূর্ণসংখ্যা):", np.clip(np.rint(demand), 0, None).astype(np.int64))

products = np.array(["চা", "ডাল", "চাল", "তেল", "চিনি"])
print("ফেরত দিয়ে:", rng.choice(products, size=6))
print("ফেরত ছাড়া:", rng.choice(products, size=3, replace=False))

pay = rng.choice(["COD", "bKash", "card"], size=10_000, p=[0.6, 0.3, 0.1])
kinds, counts = np.unique(pay, return_counts=True)
print(dict(zip(kinds.tolist(), (counts / pay.size).round(3).tolist())))

ids = np.arange(1, 7)
print("permutation:", rng.permutation(ids), "| মূল:", ids)
rng.shuffle(ids)                                   # জায়গাতেই বদলায়
print("shuffle-এর পরে মূল:", ids)

child_rngs = np.random.default_rng(7).spawn(2)     # আলাদা stream
print("দুই সন্তান stream:", [r.integers(0, 100, size=3).tolist() for r in child_rngs])

একই seed: [1 5 4 3 3] [1 5 4 3 3]
পরের ৫টি: [6 1 5 2 1]
random(): [0.179 0.64  0.467]
normal: [48.4 43.8 38.2 38.1 41.8 38.4 38.6]
দৈনিক অর্ডার (পূর্ণসংখ্যা): [48 44 38 38 42 38 39]
ফেরত দিয়ে: ['চা' 'চিনি' 'তেল' 'চিনি' 'ডাল' 'তেল']
ফেরত ছাড়া: ['ডাল' 'চা' 'তেল']
{'COD': 0.594, 'bKash': 0.298, 'card': 0.108}
permutation: [4 2 1 3 6 5] | মূল: [1 2 3 4 5 6]
shuffle-এর পরে মূল: [2 3 6 5 1 4]
দুই সন্তান stream: [[30, 79, 20], [68, 48, 39]]

৭ · পুনরুৎপাদন আর আলাদা stream — seed, spawn আর SeedSequence

Reproducibility-এর জন্য seed দরকার, কিন্তু যথেষ্ট নয়। NumPy-র docs-এর “compatibility policy” পাতা নিজেই শর্তগুলো লিখেছে: একই BitGenerator, একই seed, একই ক্রমে একই argument দিয়ে একই method ডাকা, একই NumPy build, একই পরিবেশ, একই যন্ত্র হলে একই ধারা পাওয়ার কথা, আর docs মেনে নেয় শর্তগুলো “খুব কড়া”। কয়েকটি ব্যবহারিক ফল:

• ডাকার আকারও ধারার অংশ। docs অনুযায়ী rng.random() পাঁচবার ডাকা আর rng.random(5) একবার ডাকা একই সংখ্যা দেবে, এমন নিশ্চয়তা নেই।
• Version বদলালে Generator-এর ধারা বদলাতে পারে। নতুন feature বা গতির জন্য Generator-এর ধারা ভাঙার অনুমতি আছে, তবে সাবধানে আর শুধু feature release-এ; ভুল সংশোধনে bugfix release-এও। উল্টো দিকে পুরোনো RandomState-এর ধারা স্থির রাখা হয়, ভুল থাকলেও আর বদলানো হবে না। তাই seed-এর পাশে NumPy-র version লিখে রাখুন, lockfile দিয়ে (পর্ব ১২)।
• মাঝখানে একটি নতুন ডাক পুরো ধারা সরিয়ে দেয়। Simulation-এর মাঝে একটি debug লাইনে rng.random() ঢোকালেই পরের সব ফল বদলে যায়। প্রতিটি আলাদা উদ্দেশ্যের জন্য আলাদা stream এর সমাধান।

আলাদা stream — কেন seed + 1 নয়

ধরুন বাজারঘরের simulation-এ দুটি অংশ: চাহিদা আর কুরিয়ারের দেরি। আপনি চান চাহিদার অংশে কোড বদলালেও কুরিয়ারের সংখ্যা না বদলাক। আবার একই simulation চার CPU core-এ ভাগ করে চালালে প্রতিটি core-এর আলাদা ধারা লাগবে। সহজ চিন্তা: core ০-এর seed 42, core ১-এর 43, এভাবে। NumPy-র parallel random পাতা এই পদ্ধতিকে স্পষ্ট ভাষায় “UNSAFE” বলে: এক run-এ ধারাগুলো আলাদা হলেও, আরেক দিন root seed 43 দিয়ে চালালে আগের run-এর core ১-এর ধারাই আবার আসে, ফলে বিভিন্ন run-এর ফল গোপনে একে অন্যের ওপর নির্ভর করে, আর গড়ে পক্ষপাত ঢোকে।

নিরাপদ পথ দুটি। প্রথমটি spawn: np.random.default_rng(7).spawn(2) একই উৎস থেকে দুটি সন্তান Generator দেয় (ওপরের কোডের শেষ লাইন)। প্রতিটি সন্তানের SeedSequence-এ মূল seed-এর সঙ্গে সন্তানের অবস্থানও মিশে থাকে, তাই docs-এর ভাষায় “অত্যন্ত উচ্চ সম্ভাবনায়” স্বাধীন ধারা আসে। Generator.spawn এসেছে তুলনামূলক নতুন NumPy-তে (রিলিজ নোট অনুযায়ী 1.25); পুরোনো version-এ একই কাজ হয় np.random.SeedSequence(7).spawn(2) দিয়ে, তারপর প্রতিটি থেকে default_rng(child). দ্বিতীয় পথ: seed হিসেবে একটি list দেওয়া, যেমন default_rng([worker_id, root_seed]); docs এটিকে seed + i-এর নিরাপদ বিকল্প হিসেবে দেখায়।

root = np.random.SeedSequence(2026)
demand_ss, courier_ss = root.spawn(2)          # প্রতিটি অংশের নিজস্ব উৎস
rng_demand = np.random.default_rng(demand_ss)
rng_courier = np.random.default_rng(courier_ss)
print(root.entropy)                            # 2026 — log-এ লিখে রাখুন

এখানে seed না দিলে SeedSequence() অপারেটিং সিস্টেম থেকে একটি বড় এলোমেলো সংখ্যা নেয়, আর root.entropy সেটি দেখায়। প্রথমবার seed ছাড়া চালিয়ে সেই সংখ্যাটি log-এ লিখে রাখলে পরে হুবহু একই run আবার পাওয়া যায়; প্রতিটি নতুন পরীক্ষায় “৪২”-ই না বসানোর একটি পরিচ্ছন্ন উপায়।

৮ · Simulation আর Monte Carlo — দেরিতে ডেলিভারির সম্ভাবনা

বাজারঘরের ম্যানেজারের প্রশ্ন: “একটি অর্ডার ক্রেতার হাতে পৌঁছাতে ৩ দিনের বেশি লাগার সম্ভাবনা কত?” ডেলিভারির সময় দুটি অংশের যোগফল। গুদামে প্যাক হতে লাগে ০, ১ বা ২ দিন, সম্ভাবনা যথাক্রমে ০.৫, ০.৩, ০.২। কুরিয়ারে লাগে ১, ২ বা ৩ দিন, সম্ভাবনা ০.৬, ০.৩, ০.১। ধরে নিচ্ছি দুটি অংশ একে অন্যের ওপর নির্ভর করে না (independent; এর সঠিক অর্থ পর্ব ৩৫-এ)।

আগে হাতে হিসাব

মোট সময় ৩ দিনের বেশি, মানে অন্তত ৪ দিন। কোন কোন জোড়ায় তা হয়?

• গুদাম ১, কুরিয়ার ৩: 0.3 × 0.1 = 0.03
• গুদাম ২, কুরিয়ার ২: 0.2 × 0.3 = 0.06
• গুদাম ২, কুরিয়ার ৩: 0.2 × 0.1 = 0.02

P(মোট > 3) = 0.03 + 0.06 + 0.02 = 0.11

স্বাধীন দুটি অংশের সম্ভাবনা গুণ হয়, আর আলাদা আলাদা জোড়া (যারা একসঙ্গে ঘটতে পারে না) যোগ হয়। এখানে হাতে হিসাব সহজ, কারণ মাত্র ৯টি জোড়া। কিন্তু অংশ অনেক হলে, বা অংশগুলো একে অন্যের ওপর নির্ভর করলে (বৃষ্টিতে কুরিয়ার ধীর), সূত্র না খুঁজে অনেকবার ঘটনাটি ঘটিয়ে দেখা সহজ। এর নাম Monte Carlo পদ্ধতি: random নমুনা দিয়ে একটি সংখ্যার অনুমান। এখানে সঠিক উত্তর 0.11 জানা আছে বলেই আমরা অনুমানের ভুল মাপতে পারব।

একটি trial (একবারের পরীক্ষা) মানে একটি কাল্পনিক অর্ডার: ওজনসহ choice দিয়ে গুদামের দিন, আরেকটি দিয়ে কুরিয়ারের দিন, যোগ, তারপর “৩-এর বেশি?” প্রশ্নের True বা False. অনেক trial-এর মধ্যে True-র ভাগই সম্ভাবনার অনুমান। নিচে seed 0 দিয়ে প্রথম ছয়টি trial একটি একটি করে দেখুন।

1rng = np.random.default_rng(0)2store = rng.choice([0, 1, 2], size=6, p=[.5, .3, .2])3courier = rng.choice([1, 2, 3], size=6, p=[.6, .3, .1])4late = store + courier > 35late.cumsum() / np.arange(1, 7)

নমুনা বাড়লে ভুল কতটা কমে

True/False-এর গড়কে অনুমান ধরলে তার ভুলের সাধারণ মাপ হলো standard error:

SE = √( p (1 − p) / n )

এখানে p আসল সম্ভাবনা (0.11) আর n trial-এর সংখ্যা। n = 100 হলে SE = √(0.11 × 0.89 / 100) = √0.000979 ≈ 0.031; n = 10,000 হলে ≈ 0.0031; n = 10,00,000 (দশ লাখ) হলে ≈ 0.0003. নিয়মটি হলো ভুল কমে 1/√n হারে: নমুনা ১০০ গুণ করলে ভুল ১০ ভাগের ১ ভাগ, আর ভুল অর্ধেক করতে নমুনা লাগে ৪ গুণ। অনেক নমুনার গড় যে আসল মানের দিকে যায়, তার পেছনের নিয়মটির নাম law of large numbers; তত্ত্বটি পর্ব ৩৬–৩৭-এ। নিচের কোডে তিনটি আকারে অনুমান চালান, আর শেষে n = 100-এর পরীক্ষাটি ২০০ বার করে দেখুন অনুমান কতটা ছড়ায়।

import numpy as np

EXACT = 0.11                      # হাতে হিসাব (অংশ ৮)
rng = np.random.default_rng(42)

def estimate(n):
    store = rng.choice([0, 1, 2], size=n, p=[0.5, 0.3, 0.2])   # গুদাম
    courier = rng.choice([1, 2, 3], size=n, p=[0.6, 0.3, 0.1]) # কুরিয়ার
    late = (store + courier) > 3                                # ৩ দিনের বেশি
    return late.mean()            # True-র ভাগ = সম্ভাবনার অনুমান

print("n │ অনুমান │ আসল থেকে দূরত্ব │ প্রত্যাশিত ভুল (SE)")
for n in (100, 10_000, 1_000_000):
    p_hat = estimate(n)
    se = np.sqrt(EXACT * (1 - EXACT) / n)                     # standard error
    print(f"{n:>9,} {p_hat:>8.4f} {abs(p_hat - EXACT):>8.4f} {se:>8.4f}")

# একই n=100 পরীক্ষা ২০০ বার: অনুমান কতটা ছড়ায়?
many = np.array([estimate(100) for _ in range(200)])
print("n=100, ২০০ বার: সর্বনিম্ন", many.min(), "সর্বোচ্চ", many.max(),
      "std", many.std().round(4))

n │ অনুমান │ আসল থেকে দূরত্ব │ প্রত্যাশিত ভুল (SE)
      100   0.0600   0.0500   0.0313
   10,000   0.1109   0.0009   0.0031
1,000,000   0.1098   0.0002   0.0003
n=100, ২০০ বার: সর্বনিম্ন 0.03 সর্বোচ্চ 0.18 std 0.0304

একই ছাঁচে মজুত ফুরিয়ে যাওয়ার সম্ভাবনাও অনুমান করা যায় (অনুশীলন ৪)। তবে মডেলটি ভুল হলে (ধরুন বৃষ্টির দিনে কুরিয়ার ধীর হয়, অথচ আমরা দুটি অংশকে স্বাধীন ধরেছি) কোটি trial-ও ভুল উত্তরকেই নিখুঁতভাবে অনুমান করবে। Monte Carlo শুধু হিসাবের ভুল কমায়, অনুমানের (assumption) ভুল নয়। Randomized algorithm হিসেবে Monte Carlo আর Las Vegas-এর পার্থক্য আছে ABCL TECH-এর DSA কোর্সের randomized algorithm পর্বে।

৯ · AI-এর যুগে — AI সহকারীর লেখা NumPy আর simulation যাচাই

“বিক্রির টেবিল থেকে প্রতিটি জেলার ভাগ বের করো” বা “ডেলিভারির দেরি simulate করো” লিখলে একটি AI সহকারী কয়েক সেকেন্ডে কোড দেবে। এই পর্বের বিষয়গুলোতেই এমন কোড সবচেয়ে বেশি চুপচাপ ভুল করে। AI সহকারীরা প্রায়ই যা করে:

• পুরোনো API লেখে: np.random.seed আর np.random.randint, NumPy 1.24-এ সরানো np.float বা np.int, বা pandas-এর সরিয়ে দেওয়া append; কখনো এমন parameter বানায় যা নেই (যেমন np.mean(..., skipna=True), যা pandas-এর, NumPy-র নয়)।
• axis উল্টো দেয়, আর বর্গাকার টেবিলে সেই ভুল কোনো error ছাড়াই পার পেয়ে যায়; keepdims ছাড়া ভাগ করে।
• x == np.nan দিয়ে ফাঁকা খোঁজে, mean-এর nan-কে fillna(0) দিয়ে “ঠিক” করে (ফাঁকাকে শূন্য বিক্রি বানিয়ে), আর warnings.filterwarnings("ignore") দিয়ে সব সতর্কবার্তা চাপা দেয়।
• NumPy আর pandas-এর std মিশিয়ে ফেলে; অল্প trial-এর Monte Carlo ফলকে নিশ্চিত সত্য হিসেবে লেখে, ভুলের মাপ ছাড়াই; seed বদলে বদলে “ভালো” ফলটি দেখায়।
• মডেলের কাজে ডেটা ভাগের আগেই পুরো ডেটার গড় বা std দিয়ে scale করে, ফলে test ডেটার তথ্য train-এ ঢুকে যায় (data leakage; মডেলের পর্বগুলোতে বিস্তারিত), আর train ডেটাতেই মাপা ফল দেখায়।

AI-এর কোড হাতে পেলে — যাচাইয়ের চেকলিস্ট

• প্রতিটি aggregation-এর পর shape ছেপে দেখেছি; ফলের shape থেকে axis-এর ঘরটি সত্যিই বাদ গেছে
• একটি ছোট টেবিলে ২–৩টি সংখ্যা হাতে হিসাব করে মিলিয়েছি, আর একটি অ-বর্গাকার টেবিলে পরীক্ষা করেছি
• ফাঁকা ঘর np.isnan দিয়ে গুনেছি; nan বাদ দেওয়া বা শূন্য বসানো, কোনটি ঠিক, তা ডেটার অর্থ দেখে ঠিক করেছি
• কোনো global filterwarnings("ignore") নেই; সতর্কবার্তা সামলানো হয়েছে শুধু নির্দিষ্ট লাইনে
• random-এ default_rng(seed), seed আর NumPy version লেখা আছে; কয়েকটি seed-এ চালিয়ে ফলের ছড়ানোও দেখেছি
• Monte Carlo-র ফলের সঙ্গে n আর standard error আছে; নিজের version-এর অফিসিয়াল docs-এ প্রতিটি অচেনা function মিলিয়েছি; আসল ক্রেতার ডেটা AI-তে পেস্ট করিনি