সমাধান — অধ্যায় ৮.২ · Simulation Study¶
অধ্যায় ফাইল:
part-8-capstone/08-02-simulation-study.md(§৭ অনুশীলনী)। সব সিমুলেশন seednp.random.default_rng(20260619)-এ চালানো; সংখ্যাগত উত্তর reproducible। প্রয়োজনীয় import:import numpy as np,from scipy import stats। মনে রাখবেন সংখ্যা-টানার ক্রম ওsizeঅভিন্ন না-হলে canonical মান হুবহু নাও মিলতে পারে (draw-order dependence, ← §৪.১)।মূল সংজ্ঞা ও নীতি। Monte-Carlo simulation — একটি পরিচিত তাত্ত্বিক ফল পুনরুৎপাদনের জন্য নকশা-করা random পরীক্ষা, যা একই সঙ্গে তত্ত্ব ও কোড উভয় verify করে (← §৪.৬)। Monte-Carlo standard error (MC SE): গড়-অনুমানে \(\approx s/\sqrt R\), অনুপাত-অনুমানে \(\approx\sqrt{\hat p(1-\hat p)/R}\); অর্ধেক করতে \(R\) চারগুণ। KS-দূরত্ব সর্বত্র \(D_n=\sup_x\lvert F_n(x)-\Phi(x)\rvert\)। common random numbers (CRN): তুলনীয় শর্তে একই random ইনপুট ব্যবহার করে পার্থক্যের ভেদ কমানো (variance reduction)। canonical মান। E1 — CLT log-log slope \(-0.4998\) (Berry–Esseen পূর্বানুমান \(-0.5\)), \(C=0.1326\), KS \(n=5\to0.0590\), \(n=40\to0.0225\), \(n=160\to0.0108\), \(n=320\to0.0071\); E2 — bootstrap coverage \(0.9130\) (MC SE \(0.0063\)), গড় CI-প্রস্থ \(0.5869\), নামমাত্র \(0.95\); E3 — optimal degree \(5\), ন্যূনতম total test error \(0.5784\), irreducible \(\sigma^2=0.4900\); E4 — consistency RMSE log-log slope \(-0.5570\) (পূর্বানুমান \(-0.5\)), asymptotic sd \(=\lambda=2.0000\), standardized-KS \(0.0180\) (\(n=200\), \(R=40000\)), bias ধনাত্মক কিন্তু \(\to0\)।
ক · ধারণাগত (conceptual)¶
সমাধান ১ (★)¶
(ক) কেন পুনরুৎপাদন একই সঙ্গে দুটি জিনিস প্রমাণ করে। একটি পরিচিত তাত্ত্বিক ফল পুনরুৎপাদন করলে একই কাজে দুটো ভিন্ন দাবি একসঙ্গে verify হয়। (i) তত্ত্বটি সংখ্যায় সত্য: যেমন E1-এ পর্যবেক্ষিত log-log slope \(-0.4998\) Berry–Esseen-ঘোষিত তাত্ত্বিক \(-0.5\)-এর সঙ্গে মেলা নিশ্চিত করে CLT-convergence প্রকৃতই \(1/\sqrt n\) হারে ঘটে। (ii) সিমুলেশন-pipeline ঠিকঠাক: যেহেতু তত্ত্বটি আগে থেকেই প্রমাণিত (← Part VII), মিল পাওয়া মানে কোড, sampling ও পরিসংখ্যান-গণনার পুরো শৃঙ্খল নির্ভুলভাবে কাজ করছে — নয়তো জানা উত্তরটা বেরোত না (← §৪.৬)।
(খ) slope \(-0.25\) বেরোলে সন্দেহ কোথায়। সন্দেহ প্রথমে কোডে, তত্ত্বে (Berry–Esseen) নয় — কারণ Berry–Esseen সুপ্রতিষ্ঠিত ও প্রমাণিত, তাই inconsistency-এর সবচেয়ে সম্ভাব্য উৎস কোডের একটি বাগ (যেমন মানক-করণে \(\sqrt n\) ভুল জায়গায়, বা KS-গণনায় ভুল রেফারেন্স-distribution, বা draw-order-জনিত সংখ্যা-দূষণ)।
(গ) কেন পরিচিত ফল বাছা হয়। একটি জানা উত্তরই একমাত্র যা তত্ত্ব ও কোড উভয়কে একসঙ্গে verify করার সুযোগ দেয় — এই built-in verification পেরিয়ে তবেই সেই সরঞ্জাম কোনো নতুন, অজানা প্রশ্নে বিশ্বাসযোগ্যভাবে প্রয়োগ করা যায়।
এক বাক্যে: একটা জানা ফল পুনরুৎপাদন একই সঙ্গে তত্ত্বকে সংখ্যায় নিশ্চিত করে ও pipeline verify করে — মিল না-এলে সন্দেহটা প্রথমে কোডের বাগে, প্রমাণিত তত্ত্বে নয়।
সমাধান ২ (★★)¶
(ক) Berry–Esseen theorem। iid \(X_1,\dots,X_n\), mean \(\mu\), variance \(\sigma^2\in(0,\infty)\), সসীম তৃতীয় absolute moment \(\rho=\mathbb E\lvert X-\mu\rvert^3<\infty\) হলে, মানক যোগফলের empirical CDF \(F_n\) ও \(N(0,1)\)-CDF \(\Phi\)-এর মধ্যে Kolmogorov দূরত্ব একটি সর্বজনীন ধ্রুবক \(C\) দিয়ে আবদ্ধ:
অর্থাৎ CLT-এর convergence কেবল ঘটে না, তার হার-ও \(1/\sqrt n\)-এ বাঁধা।
(খ) log-log slope \(-\tfrac12\)। সীমার ডান পাশে একমাত্র \(n\)-নির্ভর অংশ \(1/\sqrt n\); বাকি সব (\(C,\rho,\sigma\)) \(n\)-নিরপেক্ষ ধ্রুবক। তাই বাস্তবে KS একটি ধ্রুবক \(C'=C\rho/\sigma^3\) গুণ \(n^{-1/2}\) হারে চলে:
এটি \(\log(\mathrm{KS})\)-বনাম-\(\log n\)-এ একটি সরলরেখা, যার slope ঠিক \(-\tfrac12\) — যার সঙ্গে canonical \(-0.4998\) মেলে।
(গ) skewness কেবল ধ্রুবক বদলায়, slope নয়। Exp(1)-এর skewness \(2\) ধনাত্মক ও বড়, তাই \(\rho/\sigma^3\) বড়, ফলে \(C'\) (অর্থাৎ intercept, canonical \(C=0.1326\)) বাড়ে; কিন্তু \(n\)-এর ঘাত \(-\tfrac12\) অপরিবর্তিত থাকে — হার distribution-free, কেবল রেখার উচ্চতা distribution-নির্ভর।
এক বাক্যে: Berry–Esseen-এর একমাত্র \(n\)-নির্ভরতা \(1/\sqrt n\), তাই \(\log(\mathrm{KS})\)-বনাম-\(\log n\) slope \(-\tfrac12\) (canonical \(-0.4998\)); base-distribution-এর skewness কেবল intercept \(C\) বদলায়, slope নয়।
সমাধান ৩ (★★)¶
(ক) কেন right-skewed ডেটায় percentile bootstrap under-cover করে। base distribution Exp(1) দৃঢ়ভাবে right-skewed, তাই \(n=40\)-এর মতো মাঝারি নমুনায় নমুনা-গড়ের distribution এখনো পুরোপুরি symmetric নয় — এতে সামান্য ডান-লেজ থাকে। percentile bootstrap সরাসরি resample-distribution-এর \(2.5\) ও \(97.5\) শতাংশক কেটে CI বানায়, কিন্তু সেই পদ্ধতি এই অবশিষ্ট asymmetry ঠিকভাবে ধরতে পারে না — দুই প্রান্তকে সমান-ভাবে ঢাকে, ফলে প্রকৃত সত্য-মান কিছুটা বেশি বার বাইরে পড়ে ⇒ empirical coverage \(0.9130\), নামমাত্র \(0.95\)-এর সামান্য নিচে (← 4.9)।
(খ) \(n\) বাড়ালে। \(n\to\infty\)-এ CLT-এর ফলে নমুনা-গড় ক্রমে normal (তাই symmetric) হয়ে যায়, asymmetry মুছে যায়, তাই coverage \(0.95\)-এর দিকে উঠে আসার কথা।
(গ) কেন under-coverage প্রকৃত প্রভাব। পার্থক্য \(0.95-0.9130=0.0370\) হলো MC SE \(0.0063\)-এর প্রায় \(0.0370/0.0063\approx5.9\) গুণ — এত বড় deviation নিছক সিমুলেশন-কোলাহলে (random ওঠানামায়) ব্যাখ্যাযোগ্য নয়, তাই under-coverage-টি বাস্তব প্রভাব (← §৪.৩)।
এক বাক্যে: Exp(1)-এর skewness \(n=40\)-এ নমুনা-গড়ের distribution-কে asymmetric রাখে বলে percentile bootstrap সামান্য under-cover করে (\(0.9130\)); \(n\) বাড়লে coverage \(0.95\)-এ উঠবে, আর \(0.0370\approx5.9\,\)MC SE বলে এই ঘাটতি প্রকৃত, কোলাহল নয়।
সমাধান ৪ (★★)¶
(ক) bias–variance বিয়োজন। একটি স্থির test-বিন্দু \(x\)-এ expected prediction error তিন অংশে ভাঙে:
যেখানে \(\sigma^2\) হলো noise-জনিত irreducible ত্রুটি — যেকোনো মডেল দিয়েও কমানো অসম্ভব (এখানে noise sd \(0.7\), তাই \(\sigma^2=0.4900\))।
(খ) কেন U-আকৃতি। degree বাড়লে মডেল নমনীয়তর হয়, তাই signal-এর বাঁক ভালো ধরে ⇒ \(\operatorname{Bias}^2\) কমে; কিন্তু একই নমনীয়তা noise-ও ধরে ফেলে (over-fit) ⇒ \(\operatorname{Var}\) বাড়ে। সরল মডেল signal ধরতে পারে না (under-fit, বড় bias), জটিল মডেল noise ধরে ফেলে (over-fit, বড় variance); দুইয়ের যোগফল তাই একটি optimal degree-এ ন্যূনতম — canonical সারণিতে degree \(5\)-এ total test error \(0.5784\) (← 6.1)। Variance monotone বাড়ে, Bias² monotone কমে/সমতল হয়, তাই Total U-আকৃতির।
(গ) কেন \(\operatorname{Bias}^2\) হঠাৎ \(2\to3\)-এ পড়ে। signal \(f(x)=\sin(1.5x)+0.5x\)-এর \(\sin\)-অংশের Taylor expansion-এ কেবল বিজোড় (odd) ঘাত থাকে (\(x,x^3,x^5,\dots\)); তাই একটি জোড়-সর্বোচ্চ-ঘাতের বহুপদী নতুন কিছু ধরে না, কিন্তু বিজোড় ঘাত যোগ হওয়ামাত্র (degree \(2\to3\)) প্রথম বড় bias-হ্রাস আসে (\(0.4812\to0.0741\)), এবং ফের \(4\to5\)-এ (\(0.0742\to0.0020\))।
এক বাক্যে: \(\mathbb E[(y-\hat f)^2]=\operatorname{Bias}^2+\operatorname{Var}+\sigma^2\); degree বাড়লে Bias² কমে কিন্তু Var বাড়ে, তাই Total U-আকৃতির (min degree \(5\), \(0.5784\)), আর \(\sin\)-এর বিজোড়-ঘাত-প্রকৃতির জন্য Bias² \(2\to3\)-এ হঠাৎ পড়ে।
খ · গণনামূলক (computational)¶
সমাধান ৫ (★)¶
(ক) MC SE গণনা। coverage একটি অনুপাত, তাই প্রতিটি ডেটাসেট একটি Bernoulli ফল ("CI সত্য মান ঢাকল কি না"); অনুপাত-অনুমানের MC SE সূত্র প্রযোজ্য। \(\hat p=0.9130\), \(D=2000\):
যা canonical \(0.0063\)-এর সঙ্গে হুবহু মেলে।
(খ) ৯৫% Monte-Carlo interval। আনুমানিক interval \(\hat p\pm1.96\,\operatorname{MC\,SE}\):
নামমাত্র \(0.95\) এই interval-এর বাইরে (\(0.9253<0.95\)) — অর্থাৎ পর্যবেক্ষিত under-coverage সিমুলেশন-কোলাহলে ব্যাখ্যাযোগ্য নয়, তাৎপর্যপূর্ণ।
(গ) MC SE অর্ধেক করতে। যেহেতু \(\operatorname{MC\,SE}\propto1/\sqrt D\), MC SE অর্ধেক (\(0.00315\)) করতে \(D\)-কে চারগুণ (\(8000\)) করতে হবে, কারণ \(\sqrt4=2\) — নির্ভুলতা ব্যয়বহুল (← §৪.২)।
এক বাক্যে: \(\sqrt{0.9130\times0.0870/2000}=0.0063\) (canonical); \(95\%\) MC interval \([0.9007,0.9253]\) নামমাত্র \(0.95\)-কে বাদ দেয় ⇒ তাৎপর্যপূর্ণ under-coverage; MC SE অর্ধেক করতে \(D\) চারগুণ।
সমাধান ৬ (★★)¶
(ক) asymptotic variance \(=\lambda^2\)। MLE-এর সাধারণ asymptotic-normality ফল বলে, regularity-শর্ত পূরণে \(\sqrt n(\hat\lambda-\lambda)\xrightarrow{d}N\bigl(0,1/I(\lambda)\bigr)\), যেখানে \(I(\lambda)\) একক-পর্যবেক্ষণের Fisher information (← 4.3)। Exp(rate \(\lambda\))-এ \(I(\lambda)=1/\lambda^2\), তাই
(খ) সংখ্যাগত মান। \(\lambda=2\) বসিয়ে asymptotic sd \(=\lambda=2.0000\), যা canonical \(2.0000\)-এর সঙ্গে মেলে।
(গ) কেন \(S\)-এর sd \(\approx1\) তত্ত্ব সমর্থন করে। \(S=\sqrt n(\hat\lambda-\lambda)/\lambda\)-কে \(\lambda\) দিয়ে ভাগ করায় তার তাত্ত্বিক asymptotic sd হয় \(\lambda/\lambda=1\); \(n=200\)-এ পর্যবেক্ষিত নমুনা-sd \(1.0091\) প্রায় \(1\), তাই \(N(0,1)\)-সীমা ঠিকঠাক (সঙ্গে KS \(0.0180\)-ও ছোট)।
এক বাক্যে: Exp-rate MLE-এর asymptotic variance \(=1/I(\lambda)=\lambda^2\), তাই asymptotic sd \(=\lambda=2.0000\); standardized \(S\)-এর sd \(1.0091\approx1\) সেই \(N(0,1)\)-সীমা নিশ্চিত করে।
সমাধান ৭ (★)¶
(ক) \(n=40\to160\) স্কেলিং। Berry–Esseen বলে \(\mathrm{KS}\approx C'\,n^{-1/2}\)। \(n\) চারগুণ হলে (\(40\to160\)) \(n^{-1/2}\) ঠিক অর্ধেক হয় (\(\sqrt{40/160}=\sqrt{1/4}=1/2\)), তাই
(খ) canonical-এর সঙ্গে তুলনা। পূর্বানুমান \(0.01125\) বনাম canonical প্রকৃত \(0.0108\) — খুব কাছাকাছি। হুবহু সমান না-হওয়ার কারণ, Berry–Esseen কেবল প্রধান পদ \(n^{-1/2}\) দেয়, কিন্তু সসীম \(n\)-এ উচ্চতর-ক্রম সংশোধন (\(o(n^{-1/2})\)-পদ) থাকে, যা মাঝারি \(n\)-এ ক্ষীণ কিন্তু শূন্য নয়; তাই সামান্য ব্যবধান।
(গ) কেন সাত-বিন্দু ধারা প্রায়-সরলরেখা। \(n=5\to0.0590\) থেকে \(n=320\to0.0071\) — এই পুরো ধারায় প্রধান পদ \(n^{-1/2}\)-ই আধিপত্য করে, তাই \(\log(\mathrm{KS})\) ও \(\log n\)-এর সম্পর্ক প্রায়-linear (slope \(-0.4998\))।
এক বাক্যে: \(\mathrm{KS}\propto n^{-1/2}\)-এ \(n:40\to160\) চারগুণ হলে KS অর্ধেক (\(0.0225\to0.01125\)), যা canonical \(0.0108\)-এর কাছাকাছি — সামান্য ব্যবধান \(o(n^{-1/2})\)-সংশোধনের জন্য, আর প্রধান পদের আধিপত্যে সাত-বিন্দু log-log ধারা প্রায়-সরলরেখা।
গ · প্রমাণ/ডিজাইন-ভিত্তিক (proof / design-based)¶
সমাধান ৮ (★★)¶
দাবি ও নকশা। তত্ত্ব: iid \(X_i\sim N(\mu,\sigma^2)\) হলে \(\dfrac{(n-1)s^2}{\sigma^2}\sim\chi^2_{n-1}\), যেখানে \(s^2\) নমুনা-ভেদ (← 4.1)। এই দাবি verification-এর Monte-Carlo পরীক্ষা:
(ক) নকশা।
- স্থির factor: \(\mu=0\), \(\sigma=1\), \(n=10\) (একটিমাত্র স্থির নমুনা-আকার; এখানে কোনো factor বদলানো হয় না, কেবল একটি distribution-দাবি verify)।
- replication-সংখ্যা: \(R=50000\) বড় রাখা হয়, যাতে empirical distribution-এর MC-কোলাহল (\(\sim1/\sqrt R\)) ছোট থাকে ও KS-অনুমান নির্ভরযোগ্য হয়।
- প্রতি replication-এ গণনা: একটি নমুনা \(X_1,\dots,X_n\) টেনে \(s_r^2\) (ddof \(=1\)) বের করে \(T_r=(n-1)s_r^2/\sigma^2\)।
- মিলিয়ে-দেখা: সংগৃহীত \(\{T_r\}\)-এর empirical distribution-কে \(\chi^2_{n-1}\)-এর সঙ্গে দুইভাবে মেলানো যায় — (i) KS-দূরত্ব stats.kstest(T, 'chi2', args=(n-1,)).statistic ছোট কি না; এবং (ii) moment-মিল: \(\mathbb E[\chi^2_{n-1}]=n-1\) ও \(\operatorname{Var}=2(n-1)\), তাই \(\bar T\approx n-1=9\) ও নমুনা-ভেদ \(\approx2(n-1)=18\) কি না।
import numpy as np
from scipy import stats
rng = np.random.default_rng(20260619)
mu, sigma, n, R = 0.0, 1.0, 10, 50_000
X = rng.normal(mu, sigma, size=(R, n)) # iid N(0,1) নমুনা
s2 = X.var(axis=1, ddof=1) # নমুনা-ভেদ
T = (n - 1) * s2 / sigma**2 # তত্ত্বে ~ chi2_{n-1}
print("mean(T) =", round(T.mean(), 4), "(তত্ত্ব n-1 =", n-1, ")")
print("var(T) =", round(T.var(ddof=0), 4), "(তত্ত্ব 2(n-1) =", 2*(n-1), ")")
print("KS(T, chi2_{n-1}) =", round(stats.kstest(T, 'chi2', args=(n-1,)).statistic, 4))
(খ) কেন বীজ স্থির জরুরি। বীজ default_rng(20260619) স্থির না-রাখলে ফল পুনরুৎপাদনযোগ্য নয় — অন্য কেউ (বা ভিন্ন সময়ে একই কোড) হুবহু একই সংখ্যা পাবে না, তাই canonical-মান মেলানো বা independent verify অসম্ভব হয়ে পড়ে (← §৪.১)।
(গ) কেন এটিও pipeline-verification-এর দৃষ্টান্ত। \(\chi^2_{n-1}\)-ফল একটি প্রমাণিত তত্ত্ব; empirical distribution এর সঙ্গে মিললে একই কাজে তত্ত্বটি সংখ্যায় নিশ্চিত হয় এবং নমুনা-ভেদ-গণনা কোড (ddof, স্কেলিং) নির্ভুল বলে verify হয় — ঠিক "জানা ফল পুনরুৎপাদন করে pipeline verify"-এর নীতি (← §৪.৬)।
এক বাক্যে: \(\mu,\sigma,n\) স্থির রেখে বড় \(R\)-এ \(T_r=(n-1)s_r^2/\sigma^2\) সংগ্রহ করে \(\chi^2_{n-1}\)-এর সঙ্গে KS বা moment-মিলে verify করলে তত্ত্ব ও নমুনা-ভেদ-কোড উভয়ই একসঙ্গে পরীক্ষিত হয় — বীজ স্থির থাকায় ফল পুনরুৎপাদনযোগ্য।
সমাধান ৯ (★★)¶
(ক) কেন এটি variance reduction (common random numbers)। ধরা যাক দুই degree-এর total error \(A\) ও \(B\)। আগ্রহ তাদের পার্থক্য \(A-B\)-এ (কোন degree ভালো)। ভেদের নিয়ম:
দুই degree-এ একই noise realisation ব্যবহার করলে \(A\) ও \(B\) একই random ইনপুট ভাগ করে, তাই তারা ধনাত্মকভাবে সহ-সম্পর্কিত (\(\operatorname{Cov}(A,B)>0\)); ফলে \(\operatorname{Var}(A-B)\) কমে — অর্থাৎ degree-তুলনা তীক্ষ্ণতর হয়। এটিই common random numbers, একটি variance-reduction কৌশল (← §৪.৪)।
(খ) independent noise হলে কী ক্ষতি। প্রতিটি degree-এ independent noise ব্যবহার করলে \(\operatorname{Cov}(A,B)=0\), তাই পার্থক্যের ভেদ বাড়ে; U-আকৃতি বক্ররেখা অমসৃণ হয়ে যায়, আর দুই কাছাকাছি degree-এর প্রকৃত পার্থক্য অতিরিক্ত কোলাহলে ঢেকে যেতে পারে ⇒ optimal-degree সিদ্ধান্ত অস্থির।
(গ) কেন এটি factor-বিচ্ছিন্নকরণের সরাসরি প্রয়োগ। noise স্থির রাখায় দুই degree-এর total-error-এর পরিবর্তন কেবল degree-এর ফল, random-ভিন্নতার নয় — অর্থাৎ একমাত্র পরিবর্তনশীল factor থাকে degree, যা এক-factor (one-factor-at-a-time) নকশার শর্ত পূরণ করে।
এক বাক্যে: একই noise ব্যবহারে \(\operatorname{Cov}(A,B)>0\) হওয়ায় \(\operatorname{Var}(A-B)\) কমে (common random numbers), তাই degree-তুলনা তীক্ষ্ণ ও optimal-degree স্থির; independent noise হলে বক্ররেখা কোলাহলে দুলত — এটি degree-কে একমাত্র factor রাখার সরাসরি প্রয়োগ।
সমাধান ১০ (★★★)¶
দাবি। Exp-rate MLE \(\hat\lambda=1/\bar X\)-এর asymptotic ভেদ \(\lambda^2\), অর্থাৎ \(\sqrt n(\hat\lambda-\lambda)\xrightarrow{d}N(0,\lambda^2)\)।
(ক) \(\bar X\)-এর CLT। \(X_i\sim\text{Exp}(\lambda)\)-এ \(\mathbb E[X]=1/\lambda\) ও \(\operatorname{Var}(X)=1/\lambda^2\)। iid গড়ের CLT সরাসরি দেয় (← 3.4):
(খ) delta method। MLE \(\hat\lambda=g(\bar X)\) যেখানে \(g(t)=1/t\), তাই \(g'(t)=-1/t^2\); মূল্যায়ন-বিন্দু \(t=\mathbb E[X]=1/\lambda\)-তে \(g'(1/\lambda)=-1/(1/\lambda)^2=-\lambda^2\)। delta method বলে একটি মসৃণ রূপান্তরের asymptotic ভেদ \(=(g')^2\times(\text{original variance})\):
(গ) কেন \(\lambda^2=1/I(\lambda)\) ও slope \(-\tfrac12\)। Exp(rate \(\lambda\))-এর Fisher information \(I(\lambda)=1/\lambda^2\), তাই \(1/I(\lambda)=\lambda^2\) — delta-method-জাত ভেদ ঠিক Cramér–Rao নিম্নসীমার সমান, অর্থাৎ MLE asymptotically efficient। যেহেতু asymptotic sd \(=\lambda\), RMSE \(\approx\lambda/\sqrt n\propto n^{-1/2}\), তাই \(\log(\mathrm{RMSE})\)-বনাম-\(\log n\)-এ slope \(-\tfrac12\) (consistency, canonical \(-0.5570\))। \(\blacksquare\)
এক বাক্যে: \(\bar X\)-এর CLT-তে \(g(t)=1/t\)-এর delta method বসালে \((g')^2\cdot\tfrac1{\lambda^2}=\lambda^4\cdot\tfrac1{\lambda^2}=\lambda^2\) পাওয়া যায়, যা \(1/I(\lambda)\)-এর সমান (efficient), আর RMSE \(\propto n^{-1/2}\) বলে log-log slope \(-\tfrac12\)।
ঘ · কোডিং (coding)¶
সব স্নিপেট seed
np.random.default_rng(20260619)-এ চালানো; সংখ্যাগত উত্তর reproducible। import:import numpy as np,from scipy import stats। সংখ্যা-টানার ক্রম ওsizeঅভিন্ন না-হলে canonical মান হুবহু নাও মিলতে পারে (draw-order dependence, ← §৪.১)।
সমাধান ১১ (★★)¶
লক্ষ্য। iid \(\text{Exp}(1)\) (mean \(1\), sd \(1\), skewness \(2\))-এর মানক নমুনা-গড় \(Z_n=\sqrt n(\bar X_n-1)\)-এর \(N(0,1)\)-এর সাপেক্ষে KS-দূরত্ব \(n\in\{5,10,20,40,80,160,320\}\)-এ (প্রতিটিতে \(R=60000\)) মেপে log-log slope বের করা — Berry–Esseen-পূর্বানুমিত \(-0.5\) পুনরুৎপাদন।
import numpy as np
from scipy import stats
rng = np.random.default_rng(20260619)
ns = np.array([5, 10, 20, 40, 80, 160, 320])
R1 = 60_000
ks = np.empty(len(ns))
for j, n in enumerate(ns):
x = rng.standard_exponential(size=(R1, n)) # iid Exp(1)
Z = np.sqrt(n) * (x.mean(axis=1) - 1.0) # মানক নমুনা-গড়
ks[j] = stats.kstest(Z, "norm").statistic # sup|F_R - Phi|
slope, intercept = np.polyfit(np.log(ns.astype(float)), np.log(ks), 1)
for n, k in zip(ns, ks):
print(f" n={n:>3} | KS={k:.4f}")
print("log-log slope =", round(slope, 4), " C = e^intercept =", round(np.exp(intercept), 4))
ফল (canonical)।
| \(n\) | KS দূরত্ব | \(C/\sqrt n\) ফিট |
|---|---|---|
| \(5\) | \(0.0590\) | \(0.0593\) |
| \(10\) | \(0.0417\) | \(0.0419\) |
| \(20\) | \(0.0285\) | \(0.0297\) |
| \(40\) | \(0.0225\) | \(0.0210\) |
| \(80\) | \(0.0148\) | \(0.0148\) |
| \(160\) | \(0.0108\) | \(0.0105\) |
| \(320\) | \(0.0071\) | \(0.0074\) |
log-log slope \(=-0.4998\) (Berry–Esseen পূর্বানুমান \(-0.5\)), \(C=e^{\text{intercept}}=0.1326\)।
(খ) মিল। KS তালিকা \(n=5\to0.0590\), \(n=40\to0.0225\), \(n=160\to0.0108\), \(n=320\to0.0071\) — canonical লক্ষ্যের সঙ্গে হুবহু; slope \(-0.4998\) তাত্ত্বিক \(-0.5\)-এর অত্যন্ত কাছে।
(গ) কেন এটি তত্ত্ব ও pipeline দুটোই verify করে। slope \(-0.5\)-এর এত কাছে আসা একই সঙ্গে Berry–Esseen-এর \(1/\sqrt n\)-হার সংখ্যায় নিশ্চিত করে এবং মানক-করণ, KS-গণনা ও sampling-কোডের নির্ভুলতা প্রমাণ করে — মিল না-এলে সন্দেহ যেত প্রথমে কোডে (← §৪.৬)। (\(R=60000\) বড় রাখা হয় যাতে KS-অনুমানের কোলাহল \(1/\sqrt R\) বৃহত্তম \(n\)-এর সংকেত \(1/\sqrt n\)-এর নিচে থাকে।)
সমাধান ১২ (★★)¶
লক্ষ্য। Exp(1) থেকে সত্য গড় \(\theta=1\); \(n=40\); \(D=2000\) ডেটাসেট; প্রতিটিতে \(B=1000\) resample-এর percentile bootstrap CI নমুনা-গড়ের ওপর — empirical coverage \(0.9130\) (MC SE \(0.0063\)) ও গড় CI-প্রস্থ \(0.5869\) পুনরুৎপাদন।
import numpy as np
rng = np.random.default_rng(20260619)
n, D, B, theta = 40, 2000, 1000, 1.0
covered = np.zeros(D, dtype=bool)
widths = np.empty(D)
for d in range(D):
data = rng.standard_exponential(size=n) # একটি Exp(1) ডেটাসেট
idx = rng.integers(0, n, size=(B, n)) # B টি bootstrap resample
boot = data[idx].mean(axis=1) # গড়ের bootstrap distribution
lo, hi = np.percentile(boot, [2.5, 97.5]) # ৯৫% percentile CI
covered[d] = (lo <= theta <= hi)
widths[d] = hi - lo
coverage = covered.mean()
mc_se = np.sqrt(coverage * (1 - coverage) / D)
print("coverage =", round(coverage, 4)) # 0.9130
print("MC SE =", round(mc_se, 4)) # 0.0063
print("mean width =", round(widths.mean(), 4)) # 0.5869
ফল (canonical)।
| রাশি | সিমুলেশন | তত্ত্ব/নামমাত্র |
|---|---|---|
| empirical coverage | \(0.9130\) | \(0.95\) |
| Monte-Carlo SE | \(0.0063\) | \(\sqrt{0.9130\times0.0870/2000}\) |
| গড় CI-প্রস্থ | \(0.5869\) | — |
(গ) কেন coverage \(0.95\)-এর নিচে ও তা তাৎপর্যপূর্ণ। Exp(1) right-skewed, তাই \(n=40\)-এ নমুনা-গড়ের distribution এখনো সামান্য asymmetric, যা percentile bootstrap ঠিকভাবে ধরতে না-পারায় সামান্য under-cover (\(0.9130\))। পার্থক্য \(0.95-0.9130=0.0370\) MC SE \(0.0063\)-এর প্রায় \(5.9\) গুণ ⇒ এই under-coverage প্রকৃত প্রভাব, সিমুলেশন-কোলাহল নয় (← §৪.৩)।
সমাধান ১৩ (★★)¶
লক্ষ্য। signal \(f(x)=\sin(1.5x)+0.5x\) on \([-3,3]\), noise sd \(\sigma=0.7\) (\(\sigma^2=0.4900\)); polynomial regression degree \(d=1,\dots,12\); \(n_{\text{train}}=30\) (স্থির design); \(R=500\) প্রশিক্ষণ-সেট, প্রতিটি degree-এ একই noise (common random numbers)। প্রতি degree-এ \(\operatorname{Bias}^2\), \(\operatorname{Var}\), Total বের করে optimal degree \(5\) (Total \(0.5784\)) ও U-আকৃতি পুনরুৎপাদন।
import numpy as np
SEED = 20260619
sigma, n_train, R3 = 0.7, 30, 500
degrees = np.arange(1, 13)
x_train = np.linspace(-3, 3, n_train) # স্থির design
x_test = np.linspace(-3, 3, 200) # ঘন টেস্ট-গ্রিড
signal = lambda x: np.sin(1.5 * x) + 0.5 * x
f_train, f_test = signal(x_train), signal(x_test)
bias2, variance, total = (np.empty(len(degrees)) for _ in range(3))
for k, d in enumerate(degrees):
rng_d = np.random.default_rng(SEED) # প্রতি degree-এ একই noise (CRN)
Vd, Vt = np.vander(x_train, d + 1), np.vander(x_test, d + 1)
preds = np.empty((R3, len(x_test)))
for r in range(R3):
y = f_train + sigma * rng_d.standard_normal(n_train)
coef, *_ = np.linalg.lstsq(Vd, y, rcond=None)
preds[r] = Vt @ coef
mean_pred = preds.mean(axis=0)
bias2[k] = np.mean((mean_pred - f_test) ** 2)
variance[k] = np.mean(preds.var(axis=0))
total[k] = bias2[k] + variance[k] + sigma**2
best = int(degrees[np.argmin(total)])
for d, b, v, t in zip(degrees, bias2, variance, total):
print(f" {d:>2} | Bias^2={b:.4f} | Var={v:.4f} | Total={t:.4f}")
print("optimal degree =", best, " min total =", round(total.min(), 4))
ফল (canonical)।
| deg | Bias² | Var | Total |
|---|---|---|---|
| \(1\) | \(0.4812\) | \(0.0305\) | \(1.0016\) |
| \(2\) | \(0.4812\) | \(0.0450\) | \(1.0162\) |
| \(3\) | \(0.0741\) | \(0.0598\) | \(0.6239\) |
| \(4\) | \(0.0742\) | \(0.0737\) | \(0.6378\) |
| \(5\) | \(0.0020\) | \(0.0865\) | \(0.5784\) (min) |
| \(6\) | \(0.0019\) | \(0.1001\) | \(0.5920\) |
| \(7\) | \(0.0002\) | \(0.1131\) | \(0.6033\) |
| \(8\) | \(0.0004\) | \(0.1288\) | \(0.6192\) |
| \(9\) | \(0.0004\) | \(0.1455\) | \(0.6359\) |
| \(10\) | \(0.0004\) | \(0.1641\) | \(0.6545\) |
| \(11\) | \(0.0004\) | \(0.1826\) | \(0.6731\) |
| \(12\) | \(0.0006\) | \(0.2041\) | \(0.6947\) |
irreducible \(\sigma^2=0.4900\); optimal degree \(=5\), ন্যূনতম total test error \(=0.5784\)। Bias² monotone কমে (বিশেষত \(2\to3\) ও \(4\to5\)-এ তীব্র পতন), Var monotone বাড়ে, তাই Total U-আকৃতির — degree \(5\)-এ ন্যূনতম।
(গ) কেন একই noise optimal-degree সিদ্ধান্ত তীক্ষ্ণ করে। সব degree-এ একই \(R=500\)টি noise realisation পুনর্ব্যবহার করায় (common random numbers) পাশাপাশি degree-এর Total-এর পার্থক্য কম ভেদ পায় (\(\operatorname{Cov}>0\)), তাই U-বক্ররেখা মসৃণ ও ন্যূনতম-বিন্দু (degree \(5\)) স্থিরভাবে চিহ্নিত হয় — independent noise হলে কাছাকাছি degree-গুলোর তুলনা কোলাহলে দুলত (← §৪.৪)।
সমাধান ১৪ (★★)¶
লক্ষ্য। মডেল Exp(rate \(\lambda=2\)), MLE \(\hat\lambda=1/\bar X\), \(I(\lambda)=1/\lambda^2\) ⇒ asymptotic sd \(=\lambda=2\)। (a) Consistency: \(n\in\{10,40,160,640\}\), \(R=4000\)-এ bias ও RMSE, এবং RMSE log-log slope (\(-0.5570\))। (b) Asymptotic normality: \(n=200\), \(R=40000\)-এ standardized \(S=\sqrt n(\hat\lambda-\lambda)/\lambda\)-এর \(N(0,1)\)-KS (\(0.0180\))।
import numpy as np
from scipy import stats
rng = np.random.default_rng(20260619)
lam = 2.0
# (a) consistency
ns = np.array([10, 40, 160, 640]); R4 = 4000
bias = np.empty(len(ns)); rmse = np.empty(len(ns))
for j, n in enumerate(ns):
x = rng.exponential(scale=1.0 / lam, size=(R4, n)) # Exp(rate 2)
lam_hat = 1.0 / x.mean(axis=1) # MLE = 1/Xbar
bias[j] = lam_hat.mean() - lam
rmse[j] = np.sqrt(np.mean((lam_hat - lam) ** 2))
sl, ic = np.polyfit(np.log(ns.astype(float)), np.log(rmse), 1)
for n, b, r in zip(ns, bias, rmse):
print(f" n={n:>3} | bias={b:+.4f} | RMSE={r:.4f}")
print("RMSE log-log slope =", round(sl, 4)) # -0.5570
# (b) asymptotic normality
x = rng.exponential(scale=1.0 / lam, size=(40_000, 200))
lam_hat = 1.0 / x.mean(axis=1)
S = np.sqrt(200) * (lam_hat - lam) / lam # -> N(0,1)
print("asymptotic sd = lambda =", round(lam, 4)) # 2.0000
print("KS(S, N(0,1)) =", round(stats.kstest(S, "norm").statistic, 4)) # 0.0180
ফল (canonical)।
| \(n\) | bias | RMSE | \(C/\sqrt n\) ফিট |
|---|---|---|---|
| \(10\) | \(+0.2182\) | \(0.7965\) | — |
| \(40\) | \(+0.0484\) | \(0.3404\) | — |
| \(160\) | \(+0.0114\) | \(0.1639\) | — |
| \(640\) | \(+0.0028\) | \(0.0775\) | — |
RMSE log-log slope \(=-0.5570\) (consistency পূর্বানুমান \(-0.5\)); asymptotic sd \(=\lambda=2.0000\); standardized-KS \(=0.0180\) (ছোট ⇒ Normal-সীমা ঠিক)। bias ধনাত্মক কিন্তু \(\to0\) (\(n=10\to+0.2182\), \(n=640\to+0.0028\))।
(গ) কেন bias ধনাত্মক অথচ consistent। \(g(t)=1/t\) একটি convex ফাংশন, তাই Jensen-এর inequality-এ \(\mathbb E[1/\bar X]>1/\mathbb E[\bar X]=\lambda\) — সসীম \(n\)-এ \(\hat\lambda=1/\bar X\)-এর bias ধনাত্মক। কিন্তু \(\bar X\xrightarrow{p}1/\lambda\) হওয়ায় \(\hat\lambda\xrightarrow{p}\lambda\), তাই MLE consistent; bias \(n\to\infty\)-এ \(0\)-এ নামে (RMSE-ও \(\propto n^{-1/2}\), slope \(-0.5570\approx-0.5\))।
এক বাক্যে: RMSE log-log slope \(-0.5570\) (\(1/\sqrt n\)-হার) ও standardized-KS \(0.0180\) (asymptotic sd \(=\lambda=2.0000\)) একসঙ্গে MLE-এর consistency ও asymptotic normality নিশ্চিত করে; \(1/\bar X\)-এর ধনাত্মক bias Jensen-জাত, কিন্তু \(\bar X\to1/\lambda\) ব