সমাধান — অধ্যায় ৬.৭ · Density Estimation, EM & Gaussian Mixture Models¶

অধ্যায় ফাইল: part-6-statistical-ml/06-07-density-em-gmm.md (§৭ অনুশীলনী)। চলমান dataset — seed-সূচক 20260619: একটি 3-component 2D Gaussian mixture, $n=600$, প্রকৃত weights $[0.40,0.35,0.25]$, অসম/তির্যক (unequal/tilted) covariance। $70/30$ split প্রয়োজন নেই (unsupervised — পুরো data fit হয়)। মূল notation: $p(x)=\sum_{k=1}^K\pi_k\,\mathcal N(x;\mu_k,\Sigma_k)$; responsibility $\gamma_{ik}=\frac{\pi_k\,\mathcal N(x_i;\mu_k,\Sigma_k)}{\sum_l\pi_l\,\mathcal N(x_i;\mu_l,\Sigma_l)}$; M-step $\pi_k=\frac1n\sum_i\gamma_{ik}$, $\mu_k=\frac{\sum_i\gamma_{ik}x_i}{\sum_i\gamma_{ik}}$, $\Sigma_k=\frac{\sum_i\gamma_{ik}(x_i-\mu_k)(x_i-\mu_k)^\top}{\sum_i\gamma_{ik}}$; ELBO $\mathcal L(q,\theta)=\mathbb E_q[\log p(x,z\mid\theta)]-\mathbb E_q[\log q(z)]$; BIC $=-2\ell(\hat\theta)+p\log n$。 canonical সংখ্যা। GMM($3$) BIC $4828.8$, per-sample log-likelihood $-3.933$, আনুমানিক weights $[0.405,0.349,0.246]$; ARI $0.97$ ($>$ k-means $0.914$)। BIC by $K$: $1\to5626.7$, $2\to5111.7$, $3\to\mathbf{4828.8}$ (সর্বনিম্ন), $4\to4857.5$, $5\to4890.9$, $6\to4925.9$। একটি ambiguous বিন্দুর responsibility $\approx[0.864,0.002,0.133]$। §ঘ-এর runnable script (seed-সূচক 20260619) উপরের সংখ্যাগুলো পুনরুৎপাদন করে (data-recipe/সংস্করণ-ভেদে শেষ অঙ্কে সামান্য হেরফের সম্ভব; গল্প অপরিবর্তিত)। গাণিতিক প্রশ্নের (৮, ৯, ১১, ১৪) সংখ্যাগুলো স্বয়ংসম্পূর্ণ — অধ্যায়-dataset থেকে স্বাধীন। সব log natural log ($\ln$)।

ক · ধারণাগত (conceptual)¶

সমাধান ১ (★)¶

(ক) দুই-ধাপের generative গল্প; latent variable $z_i$। formula $p(x)=\sum_k\pi_k\mathcal N(x;\mu_k,\Sigma_k)$ ঠিক একটা দুই-ধাপের নমুনা-প্রক্রিয়া বর্ণনা করে: (i) প্রথমে একটা component-label টানা হয় $z_i\sim\text{Categorical}(\pi_1,\dots,\pi_K)$ — অর্থাৎ $\pi_k$ probability-তে component $k$ বাছা হয়; (ii) তারপর সেই বাছা component থেকে বিন্দুটা টানা হয় $x_i\sim\mathcal N(\mu_{z_i},\Sigma_{z_i})$। এখানে latent (লুকানো) variable হলো এই $z_i$ — কোন Gaussian component থেকে বিন্দু $x_i$ এসেছে তার সূচক।

(খ) incomplete data — কেন্দ্রীয় অসুবিধা। data-তে আমরা শুধু $x_i$ দেখি, কিন্তু $z_i$ দেখি না — এই অর্ধেক-দেখা data-কে বলে incomplete (অসম্পূর্ণ) data (বা missing-label data)। যদি $z_i$ জানা থাকত, সমস্যা তুচ্ছ হতো: প্রতিটি component-এর বিন্দুগুলো আলাদা করে নিয়ে সাধারণ Gaussian-MLE (৪.৩) চালালেই $\mu_k,\Sigma_k$ পাওয়া যেত, আর $\pi_k$ = সেই component-এ পড়া বিন্দুর ভগ্নাংশ। কিন্তু $z_i$ অজানা বলে component-assignment আর parameter পরস্পর-নির্ভর হয়ে যায় — parameter জানলে assignment বের হয় (posterior), আবার assignment জানলে parameter বের হয় (MLE)। এই মুরগি-ডিম চক্রই EM-এর প্রয়োজন তৈরি করে।

(গ) তিন parameter-এর ভৌত অর্থ। $\pi_k$ = mixing weight — component $k$-এর জনসংখ্যা-অনুপাত/prior probability (যত বড়, তত বেশি বিন্দু সেখান থেকে আসে), $\sum_k\pi_k=1$। $\mu_k$ = component $k$-এর কেন্দ্র (গড় অবস্থান)। $\Sigma_k$ = component $k$-এর আকার, অভিমুখ ও বিস্তার (covariance) — diagonal হলে অক্ষ-সমান্তরাল উপবৃত্ত, off-diagonal $\ne0$ (tilted) হলে তির্যক উপবৃত্ত; বড় eigenvalue = সেই দিকে বেশি ছড়ানো।

সমাধান ২ (★)¶

দিক	KDE (kernel density)	mixture model (GMM)
(ক) component-সংখ্যা	$n$ — প্রতি data-বিন্দুতে এক kernel; data বাড়লে সংখ্যাও বাড়ে ⇒ nonparametric	ধ্রুবক $K$ (অল্প) — data-আকার-নিরপেক্ষ ⇒ parametric
(খ) প্রতি মূল্যায়নে term	$n$টি (পুরো data সঞ্চয়/স্ক্যান করতে হয়)	মাত্র $K$টি Gaussian — বড় $n$-এ অনেক সাশ্রয়ী ও compact

(ক) parametric বনাম nonparametric। KDE-তে "component"-সংখ্যা হুবহু data-সংখ্যা $n$ (প্রতিটি বিন্দুতে একটা kernel $K_h(x-x_i)$ বসে), তাই data বাড়লে মডেলের জটিলতাও বাড়ে — এটাই nonparametric (parameter-সংখ্যা data-র সাথে বাড়ে)। mixture-এ component-সংখ্যা একটা স্থির ছোট $K$, parameter-সংখ্যা $n$-নিরপেক্ষ — তাই parametric।

(খ) মূল্যায়ন-খরচ। নতুন $x$-এ density বের করতে KDE-কে $n$টি kernel যোগ করতে হয় (সব training-বিন্দু লাগে, বড় $n$-এ ধীর ও স্মৃতি-ভারী); fitted GMM-কে মাত্র $K$টি Gaussian যোগ করলেই হয় (data আর লাগে না)। তাই বড় $n$-এ GMM মূল্যায়নে ও সঞ্চয়ে অনেক সাশ্রয়ী।

(গ) কখন mixture স্বাভাবিক। যখন data প্রকৃতপক্ষে অল্প-কয়েকটা গুচ্ছ/subpopulation থেকে আসে (চলমান উদাহরণের $3$টি Gaussian-এর মতো) — তখন mixture শুধু সাশ্রয়ীই নয়, ব্যাখ্যামূলকও: প্রতিটি component একটা অর্থপূর্ণ গোষ্ঠী (যেমন গ্রাহক-শ্রেণি, কোষ-ধরন), আর $\pi_k,\mu_k,\Sigma_k$ সেই গোষ্ঠীর আকার-কেন্দ্র-বিস্তার বলে দেয়।

সমাধান ৩ (★)¶

(ক) responsibility = Bayes posterior। $\gamma_{ik}$ হুবহু posterior probability $P(z_i{=}k\mid x_i)$ — Bayes-এর সূত্রে (২.২): $$ \gamma_{ik}=P(z_i{=}k\mid x_i)=\frac{P(z_i{=}k)\,p(x_i\mid z_i{=}k)}{p(x_i)}=\frac{\overbrace{\pi_k}^{\text{prior}}\,\overbrace{\mathcal N(x_i;\mu_k,\Sigma_k)}^{\text{likelihood}}}{\underbrace{\sum_l\pi_l\mathcal N(x_i;\mu_l,\Sigma_l)}_{\text{evidence }=\,p(x_i)}}. $$ prior = $\pi_k$ (component $k$-এর জনসংখ্যা-অনুপাত), likelihood = $\mathcal N(x_i;\mu_k,\Sigma_k)$ (component $k$ থেকে $x_i$ আসার ঘনত্ব), evidence = হর $p(x_i)$ (সব component-এর marginal ঘনত্ব)।

(খ) $\sum_k\gamma_{ik}=1$ ⇒ soft assignment। যেহেতু $\gamma_{ik}$ একটা posterior probability-বণ্টন (component-গুলোর উপর), $\sum_{k=1}^K\gamma_{ik}=1$ (হর সব লব-এর যোগফল, তাই normalize হয়ে $1$)। তাই প্রতিটি বিন্দু সব component-এ ভগ্নাংশে ভাগ হয়ে যায় (যেমন $[0.864,0.002,0.133]$) — কোনো একটিতে সম্পূর্ণ ($1$) নয়। এটাই soft assignment: বিন্দুটি "$86.4\%$ component-$1$, $13.3\%$ component-$3$"। k-means-এর hard assignment বিন্দুটিকে ঠিক একটা cluster-এ ($0/1$) দেয় — অনিশ্চয়তা ফেলে দেয়।

(গ) ambiguous বিন্দু $[0.864,0.002,0.133]$। এই বিন্দু মূলত component-$1$-এর ($86.4\%$), সামান্য component-$3$-এ ($13.3\%$), component-$2$-তে কার্যত নয় ($0.2\%$) — অর্থাৎ এটা component-$1$ ও $3$-এর সীমানার কাছে (component-$2$ থেকে দূরে)। hard assignment হলে এটা সর্বোচ্চ-$\gamma$-র component-$1$ পেত; কিন্তু "$13.3\%$ সম্ভাবনায় এটা আসলে component-$3$-এর" — এই সীমানা-অনিশ্চয়তা hard-এ হারিয়ে যেত। soft assignment এই অনিশ্চয়তা সততার সাথে ধরে রাখে।

সমাধান ৪ (★★)¶

(ক) log-sum কেন closed-form আটকায়। log-likelihood $$ \ell(\theta)=\sum_{i=1}^n\log\Big(\sum_{k=1}^K\pi_k\,\mathcal N(x_i;\mu_k,\Sigma_k)\Big). $$ একক Gaussian-এ (৪.৩) $\log\mathcal N$ একটা পরিষ্কার quadratic (কারণ $\log\exp(\cdot)=(\cdot)$ — $\log$ আর $\exp$ বাতিল হয়), তাই derivative শূন্য করলে রৈখিক সমীকরণ ⇒ closed-form। কিন্তু এখানে $\log$-এর ভেতরে যোগফল ($\log\sum_k$) থাকায় সেই বাতিল ঘটে না: $\mu_k$-সাপেক্ষে differentiate করলে $$ \frac{\partial\ell}{\partial\mu_k}=\sum_i\frac{\pi_k\,\nabla_{\mu_k}\mathcal N(x_i;\mu_k,\Sigma_k)}{\sum_l\pi_l\mathcal N(x_i;\mu_l,\Sigma_l)}=\sum_i\underbrace{\frac{\pi_k\mathcal N(x_i;\mu_k,\Sigma_k)}{\sum_l\pi_l\mathcal N_l}}{=\,\gamma(x_i-\mu_k), $$ যেখানে হরে পুরো mixture }}\,\Sigma_k^{-1$\sum_l\pi_l\mathcal N_l$ থেকে যায় — অর্থাৎ $\mu_k$-এর সমীকরণে $\gamma_{ik}$ লুকিয়ে আছে, আর $\gamma_{ik}$ নিজেই সব $\mu_l,\Sigma_l,\pi_l$-এর উপর নির্ভরশীল। তাই সমীকরণগুলো coupled ও nonlinear — কোনো একটা $\mu_k$ আলাদা করে বের করা যায় না, closed-form নেই।

(খ) complete data হলে সহজ। যদি প্রতিটি $z_i$ জানা থাকত, complete-data log-likelihood হতো $$ \ell_c(\theta)=\sum_i\log\big(\pi_{z_i}\mathcal N(x_i;\mu_{z_i},\Sigma_{z_i})\big)=\sum_i\big[\log\pi_{z_i}+\log\mathcal N(x_i;\mu_{z_i},\Sigma_{z_i})\big]. $$ এখানে $\log$-এর ভেতরে আর যোগফল নেই (একটা product-এর log, log-sum নয়) — তাই বিন্দুগুলোকে component অনুযায়ী দলবদ্ধ করে প্রতিটি component আলাদা সাধারণ Gaussian-MLE (৪.৩, closed-form) করলেই $\mu_k,\Sigma_k$ মেলে, আর $\pi_k=$ (component $k$-এ পড়া বিন্দুর ভগ্নাংশ)। log-sum আর log-product-এর পার্থক্যই কঠিন-বনাম-সহজের মূল।

(গ) EM কীভাবে পাশ কাটায়। EM অজানা $z_i$-কে প্রথমে E-step-এ soft responsibility $\gamma_{ik}$ দিয়ে "পূরণ" করে (posterior), তারপর M-step-এ ঠিক সেই $\gamma$-ওজনে সহজ weighted Gaussian-MLE করে — অর্থাৎ log-sum-এর কঠিন সমস্যাকে দুটো সহজ ধাপে (posterior + weighted-MLE) ভেঙে পুনরাবৃত্তভাবে সমাধান করে।

সমাধান ৫ (★★)¶

decomposition: $\log p(x\mid\theta)=\mathcal L(q,\theta)+\mathrm{KL}\big(q(z)\,\Vert\,p(z\mid x,\theta)\big)$, $\mathrm{KL}\ge0$।

(ক) ELBO একটা lower bound। যেহেতু সর্বদা $\mathrm{KL}(\cdot\Vert\cdot)\ge0$ (Gibbs/Jensen), পরিচয় থেকে $\log p(x\mid\theta)=\mathcal L(q,\theta)+\mathrm{KL}\ge\mathcal L(q,\theta)$ — তাই $\mathcal L$ সবসময় $\log p(x)$ (log-evidence)-এর নিচে বা সমান, অর্থাৎ একটা lower bound। নাম "evidence lower bound (ELBO)" তাই যথাযথ: $\mathcal L$ হলো evidence-এর নিচের একটা সীমা।

(খ) E-step bound-কে tight করে। $\theta$ স্থির রেখে $q(z)\leftarrow p(z\mid x,\theta)$ (= responsibility) বসালে $\mathrm{KL}\big(q\,\Vert\,p(z\mid x,\theta)\big)=\mathrm{KL}\big(p(z\mid x,\theta)\,\Vert\,p(z\mid x,\theta)\big)=0$ (একই distribution-এর মধ্যে KL শূন্য)। তাই $\log p(x\mid\theta)=\mathcal L(q,\theta)+0=\mathcal L(q,\theta)$ — bound ঠিক বর্তমান $\theta$-তে $\log p(x)$-কে স্পর্শ করে (আঁটসাঁট/tight)। অর্থাৎ E-step lower bound-কে উপরে ঠেলে current likelihood-এ ছোঁয়ায়।

(গ) M-step likelihood বাড়ায় (monotone)। M-step ($q$ স্থির রেখে $\mathcal L$-কে $\theta$-সাপেক্ষে বাড়ানো) ELBO $\mathcal L$ বাড়ায় (বা স্থির রাখে)। এখন: E-step-এর পর $\log p(x\mid\theta^{\text{old}})=\mathcal L$ (tight), আর সবসময় $\log p(x\mid\theta^{\text{new}})\ge\mathcal L(q,\theta^{\text{new}})$ (কারণ $\mathrm{KL}\ge0$)। M-step-এ $\mathcal L$ বাড়লে $\log p(x\mid\theta^{\text{new}})\ge\mathcal L(q,\theta^{\text{new}})\ge\mathcal L(q,\theta^{\text{old}})=\log p(x\mid\theta^{\text{old}})$ — তাই প্রতি iteration-এ log-likelihood কখনো কমে না (monotone↑)। মূল চালিকা: Jensen-এর অসমতা থেকে আসা $\mathrm{KL}\ge0$ (পূর্ণ শৃঙ্খল — সমাধান ১২)।

সমাধান ৬ (★★)¶

(ক) tilted-elliptical-এ k-means কেন পিছিয়ে। k-means প্রতিটি বিন্দুকে নিকটতম centroid-এ assign করে — তার সিদ্ধান্ত-সীমানা সবসময় centroid-জোড়ার লম্ব-সমদ্বিখণ্ডক (Voronoi), যা কার্যত প্রতিটি cluster-কে spherical (গোলকীয়) ও সমান-আকারের ধরে নেয় (Euclidean দূরত্ব)। চলমান উদাহরণে covariance tilted/unequal — cluster-গুলো তির্যক, লম্বাটে, ভিন্ন-আকারের উপবৃত্ত। এমন cluster-এ k-means-এর সোজা-গোল সীমানা প্রকৃত উপবৃত্তীয় সীমানা থেকে সরে যায়, তাই overlapping/তির্যক অঞ্চলের বিন্দু ভুল cluster পায় ⇒ ARI কম ($0.914$)। GMM প্রতিটি $\Sigma_k$ আলাদা শিখে উপবৃত্তীয় সীমানা টানে, তাই ভালো করে (ARI $0.97$)।

(খ) কোন দুই সরলীকরণে GMM = k-means। (i) সব $\Sigma_k=\sigma^2 I$ ধরা (সমান, গোলকীয় covariance — সব দিকে সমান বিস্তার, component-নিরপেক্ষ) এবং (ii) soft responsibility-কে hard করা (winner-take-all — $\sigma\to0$ সীমায় $\gamma_{ik}\to1$ নিকটতম-$\mu$-র জন্য, বাকিতে $0$)। এই দুই শর্তে GMM-এর E-step "নিকটতম centroid বাছা"-তে আর M-step "assign-করা বিন্দুর গড় নেওয়া"-তে পরিণত হয় — হুবহু k-means।

(গ) k-means = GMM-এর বিশেষ সীমা। GMM বেশি সাধারণ; k-means হলো GMM-এর "সমান-গোলকীয় covariance + hard assignment" বিশেষ সীমা — তাই k-means-কে "hard-spherical EM"-ও বলা হয়। GMM যা পারে (elliptical, soft, ভিন্ন-আকার), k-means তার একটা সংকীর্ণ-অনুমান উপসেট।

সমাধান ৭ (★★)¶

BIC $=-2\ell(\hat\theta)+p\log n$ (কম = ভালো)। canonical: $K{=}1\to5626.7$, $2\to5111.7$, $3\to\mathbf{4828.8}$, $4\to4857.5$, $5\to4890.9$, $6\to4925.9$।

(ক) কোন $K$, কেন। BIC যত কম তত ভালো (কম $-2\ell$ = ভালো fit, কম $p\log n$ = কম জটিলতা), তাই সর্বনিম্ন BIC-র $K$ বাছি — এখানে $K=\mathbf{3}$ ($4828.8$, সবচেয়ে ছোট মান)। (সর্বোচ্চ নয় — BIC একটা penalized-loss, minimize করার রাশি।)

(খ) "$3$-এ ভাঙা" ও $K>3$-এ BIC ফের বাড়া। প্রকৃত data $3$টি Gaussian থেকে। $K<3$-এ model আসল গঠন ধরতে পারে না — একটা বা দুটো Gaussian দিয়ে তিনটি গোষ্ঠী ব্যাখ্যা করতে গিয়ে likelihood অনেক কম (বড় $-2\ell$), তাই BIC বড় ($5626.7,5111.7$) এবং $K$ বাড়ালে দ্রুত পড়ে। $K=3$-এ model ঠিক প্রকৃত গঠন ধরে — likelihood সর্বোচ্চ-কার্যকর, BIC সর্বনিম্ন। $K>3$-এ অতিরিক্ত component প্রকৃত data-তে নতুন কিছু ব্যাখ্যা করে না — likelihood সামান্যই বাড়ে, কিন্তু প্রতিটি নতুন component অনেক parameter যোগ করে ($\pi,\mu,\Sigma$), তাই $p\log n$ penalty বেশি বাড়ে ⇒ যোগফল BIC আবার ওঠে ($4857.5,4890.9,4925.9$)। ফলে $K=3$-এ স্পষ্ট সর্বনিম্ন (V-আকার)।

(গ) দুই পদ — Occam-এর ক্ষুর। $-2\ell$ হলো misfit (fit খারাপ হলে বড়; component বাড়ালে ছোট হয়)। $p\log n$ হলো complexity penalty ($p$ = parameter-সংখ্যা, component বাড়ালে বড়; $\log n$ = নমুনা-আকার-নির্ভর কঠোরতা)। যোগফল minimize করা মানে "fit যথেষ্ট ভালো রাখো, কিন্তু যত কম parameter দিয়ে সম্ভব" — অর্থাৎ Occam-এর ক্ষুর (অপ্রয়োজনীয় জটিলতা এড়াও)-এর সংখ্যাগত রূপ। BIC-র $\log n$ factor AIC-র $2$-এর চেয়ে কঠোর, তাই BIC সাধারণত ছোট (parsimonious) model বাছে।

খ · গণনামূলক (computational)¶

সমাধান ৮ (★★)¶

১D ২-component: $\pi=(0.6,0.4)$, $\mu=(0,4)$, $\sigma=(1,1)$; $x=1.5$। ১D Gaussian $\mathcal N(x;\mu,\sigma^2)=\frac{1}{\sqrt{2\pi}\,\sigma}\exp\!\big(-\frac{(x-\mu)^2}{2\sigma^2}\big)$, এখানে $\frac{1}{\sqrt{2\pi}}=0.3989$।

(ক) component-density। $$ \mathcal N(1.5;0,1)=0.3989\,e^{-(1.5-0)^2/2}=0.3989\,e^{-1.125}=0.3989(0.32465)=\mathbf{0.12952}. $$ $$ \mathcal N(1.5;4,1)=0.3989\,e^{-(1.5-4)^2/2}=0.3989\,e^{-3.125}=0.3989(0.04394)=\mathbf{0.01753}. $$

(খ) ওজনিত-পদ ও $p(x)$। $$ \pi_1\mathcal N_1=0.6(0.12952)=\mathbf{0.07771},\qquad \pi_2\mathcal N_2=0.4(0.01753)=\mathbf{0.00701}. $$ $$ p(x)=0.07771+0.00701=\mathbf{0.08472}. $$

(গ) responsibility ও যাচাই। $$ \gamma_1=\frac{0.07771}{0.08472}=\mathbf{0.917},\qquad \gamma_2=\frac{0.00701}{0.08472}=\mathbf{0.083}. $$ $\gamma_1+\gamma_2=0.917+0.083=\mathbf{1.000}$ ✓। component-$1$ ($\mu{=}0$) এই বিন্দুর দায়িত্ব অনেক বেশি নেয় ($91.7\%$), কারণ $x{=}1.5$ তার কেন্দ্র $0$-এর অনেক কাছে (দূরত্ব $1.5$), অথচ component-$2$-এর কেন্দ্র $4$ থেকে অনেক দূরে (দূরত্ব $2.5$) — Gaussian-ঘনত্ব দূরত্বের সাথে দ্রুত ($\exp$-quadratic-এ) পড়ে বলে কাছের component বিপুলভাবে বেশি responsibility পায়।

সমাধান ৯ (★★)¶

$x=(1,2,6,8)$, $\gamma_{\cdot k}=(0.9,0.8,0.2,0.1)$, $n=4$; M-step $N_k=\sum_i\gamma_{ik}$, $\pi_k=N_k/n$, $\mu_k=\frac{\sum_i\gamma_{ik}x_i}{N_k}$।

(ক) effective count $N_k$। $$ N_k=\sum_i\gamma_{ik}=0.9+0.8+0.2+0.1=\mathbf{2.0}. $$ এটা পূর্ণসংখ্যা নয় (এখানে কাকতালীয়ভাবে $2.0$, সাধারণত ভগ্নাংশ) কারণ প্রতিটি বিন্দু এই component-এর কেবল একটা ভগ্নাংশ — $N_k$ হলো এই component-কে দায়ী-ধরা বিন্দুর "কার্যকর (effective) সংখ্যা" (soft গণনা: পূর্ণ বিন্দু নয়, responsibility-যোগফল)। এখানে $0.9+0.8+0.2+0.1=2.0$, যেন "$2$টি বিন্দুর সমতুল্য" এই component-এর।

(খ) $\pi_k$। $$ \pi_k=\frac{N_k}{n}=\frac{2.0}{4}=\mathbf{0.5}. $$

(গ) weighted mean $\mu_k$। $$ \sum_i\gamma_{ik}x_i=0.9(1)+0.8(2)+0.2(6)+0.1(8)=0.9+1.6+1.2+0.8=4.5, $$ $$ \mu_k=\frac{4.5}{2.0}=\mathbf{2.25}. $$ সাধারণ (unweighted) গড় $\frac{1+2+6+8}{4}=4.25$-এর চেয়ে এটা অনেক ছোট, কারণ weighted mean responsibility দিয়ে ওজন করে: ছোট বিন্দু $1,2$ বড় responsibility ($0.9,0.8$) পেল (তারা এই component-এর কাছের), বড় বিন্দু $6,8$ ছোট responsibility ($0.2,0.1$) পেল (তারা অন্য component-এর) — তাই গড় ছোট-বিন্দুর দিকে টানা ($2.25$), $4.25$ নয়। এটাই M-step-এর সৌন্দর্য: প্রতিটি বিন্দু তার responsibility-অনুপাতে নিজ component-এর কেন্দ্র নির্ধারণে অবদান রাখে।

সমাধান ১০ (★★)¶

ambiguous বিন্দুর responsibility $\gamma=[0.864,0.002,0.133]$।

(ক) soft assignment। তিন সংখ্যার যোগফল $0.864+0.002+0.133=0.999\approx1$ (rounding-জনিত)। এরা = বিন্দুটি প্রতিটি component থেকে আসার posterior probability $P(z{=}k\mid x)$। soft assignment হিসেবে বিন্দুটি প্রায় পুরোটা component-$1$-এর ($86.4\%$), সামান্য component-$3$-এর ($13.3\%$), component-$2$-তে কার্যত নয় ($0.2\%$) — বিন্দুটা একটা component-এ "বন্দি" নয়, বরং (মূলত দুই) component-এর মধ্যে ভগ্নাংশে বণ্টিত।

(খ) hard assignment ও হারানো তথ্য। hard (MAP/winner-take-all) assignment = সর্বোচ্চ-$\gamma$-র component, এখানে $\arg\max=$ component-$1$ ($0.864$ সর্বোচ্চ)। hard করলে যে তথ্য হারায়: "$13.3\%$ সম্ভাবনায় বিন্দুটি আসলে component-$3$-এর" — অর্থাৎ এটা component-$1$ ও $3$-এর সীমানা-সংলগ্ন/অনিশ্চিত বিন্দু, এই সীমানা-তথ্য hard label ($z=1$) সম্পূর্ণ মুছে দেয় (model আর জানে না বিন্দুটি সন্দেহজনক)।

(গ) প্রায়-সমান responsibility $[0.34,0.33,0.33]$। soft-এ এটা বলে "তিন component-এর মধ্যে প্রায়-সমান অনিশ্চয়তা" — বিন্দুটি অত্যন্ত ambiguous, কার্যত তিন cluster-এর সংযোগস্থলে। hard জোর করে component-$1$ বেছে নেবে যদিও আস্থা মাত্র $34\%$ (প্রায় random)। এমন প্রায়-সমান বিন্দুতে hard সিদ্ধান্ত বিশেষভাবে বিভ্রান্তিকর: সামান্য noise-এই সর্বোচ্চ-$\gamma$ এক component থেকে আরেকটায় লাফ দিতে পারে (অস্থির), আর "$34\%$ আস্থা"-কে "$100\%$ component-$1$" বানিয়ে মিথ্যা নিশ্চয়তা দেয়। soft assignment সততার সাথে এই অনিশ্চয়তা ধরে রাখে — downstream সিদ্ধান্তে (যেমন abstain/মানব-পর্যালোচনা) কাজে লাগে।

গ · প্রমাণভিত্তিক (proof-based)¶

সমাধান ১১ (★★★)¶

$Q(\theta)=\sum_i\sum_k\gamma_{ik}[\log\pi_k+\log\mathcal N(x_i;\mu_k,\Sigma_k)]$, $\gamma_{ik}$ স্থির; $\log\mathcal N(x_i;\mu_k,\Sigma_k)=-\tfrac12(x_i-\mu_k)^\top\Sigma_k^{-1}(x_i-\mu_k)+c$ ($c$ = $\mu_k$-নিরপেক্ষ)।

(ক) $\mu_k$-সাপেক্ষে gradient। $Q$-এ শুধু $\log\mathcal N(x_i;\mu_k,\Sigma_k)$ পদ $\mu_k$-নির্ভর। multivariate-quadratic-এর gradient (vector calculus): $\nabla_{\mu_k}\big[-\tfrac12(x_i-\mu_k)^\top\Sigma_k^{-1}(x_i-\mu_k)\big]=\Sigma_k^{-1}(x_i-\mu_k)$ ($\Sigma_k$ symmetric)। তাই $$ \frac{\partial Q}{\partial\mu_k}=\sum_{i=1}^n\gamma_{ik}\,\Sigma_k^{-1}(x_i-\mu_k)\overset{!}{=}0. $$

(খ) $\mu_k$-update সমাধান। $\Sigma_k^{-1}$ invertible (positive-definite covariance), তাই বাঁ থেকে $\Sigma_k$ গুণ করে বাদ দেওয়া যায়: $$ \sum_i\gamma_{ik}(x_i-\mu_k)=0\;\Longrightarrow\;\sum_i\gamma_{ik}x_i=\mu_k\sum_i\gamma_{ik}\;\Longrightarrow\;\boxed{\mu_k=\frac{\sum_i\gamma_{ik}x_i}{\sum_i\gamma_{ik}}}\quad\blacksquare $$ অর্থাৎ $\mu_k$ = data-র responsibility-weighted mean — প্রতিটি বিন্দু তার $\gamma_{ik}$-অনুপাতে component $k$-এর কেন্দ্র নির্ধারণে অবদান রাখে (হর $\sum_i\gamma_{ik}=N_k$ = effective count)। [একইভাবে $\frac{\partial Q}{\partial\pi_k}$ ($\sum_k\pi_k=1$ constraint-সহ, Lagrange) দেয় $\pi_k=\frac{N_k}{n}$, আর $\frac{\partial Q}{\partial\Sigma_k}$ দেয় weighted covariance — সব M-step সূত্র এই একই $Q$-maximization থেকে।]

(গ) hard-সীমা = k-means centroid। $\gamma_{ik}\in\{0,1\}$ (hard — প্রতিটি বিন্দু ঠিক একটা component-এ) হলে $\sum_i\gamma_{ik}x_i=\sum_{i:\,z_i=k}x_i$ (component $k$-এ assign-করা বিন্দুর যোগফল) আর $\sum_i\gamma_{ik}=\lvert\{i:z_i=k\}\rvert$ (সেই বিন্দুর সংখ্যা), তাই $\mu_k=$ component-$k$-এ assign-করা বিন্দুগুলোর সাধারণ গড় — হুবহু k-means-এর centroid-update (৫.৯)। এটাই প্রমাণ করে k-means = hard-EM-এর M-step (soft responsibility-কে $0/1$-এ ভেঙে)।

সমাধান ১২ (★★★)¶

decomposition: $\log p(x\mid\theta)=\mathcal L(q,\theta)+\mathrm{KL}\big(q\,\Vert\,p(z\mid x,\theta)\big)$, $\mathcal L(q,\theta)=\mathbb E_q[\log p(x,z\mid\theta)]-\mathbb E_q[\log q(z)]$, $\mathrm{KL}\ge0$ (Jensen/Gibbs)।

(ক) E-step: bound tight। iteration $t$-এ $q^{(t)}(z)=p(z\mid x,\theta^{(t)})$ নিই। তখন $$ \mathrm{KL}\big(q^{(t)}\,\Vert\,p(z\mid x,\theta^{(t)})\big)=\mathrm{KL}\big(p(z\mid x,\theta^{(t)})\,\Vert\,p(z\mid x,\theta^{(t)})\big)=0, $$ তাই decomposition থেকে $$ \mathcal L(q^{(t)},\theta^{(t)})=\log p(x\mid\theta^{(t)})-0=\log p(x\mid\theta^{(t)}).\qquad(\star) $$ অর্থাৎ E-step-এর পর lower bound ঠিক current likelihood-কে স্পর্শ করে (tight)।

(খ) M-step: ELBO বাড়ে। $\theta^{(t+1)}=\arg\max_\theta\mathcal L(q^{(t)},\theta)$ (M-step সংজ্ঞা — $q^{(t)}$ স্থির রেখে ELBO maximize)। যেহেতু $\theta^{(t+1)}$ maximizer, যেকোনো $\theta$ (বিশেষত $\theta^{(t)}$)-এর তুলনায় $$ \mathcal L(q^{(t)},\theta^{(t+1)})\ge\mathcal L(q^{(t)},\theta^{(t)}).\qquad(\star\star) $$

(গ) শৃঙ্খল জুড়ে monotonicity। যেকোনো $\theta$-তে $\mathrm{KL}\ge0$ বলে $$ \log p(x\mid\theta)=\mathcal L(q^{(t)},\theta)+\mathrm{KL}\big(q^{(t)}\,\Vert\,p(z\mid x,\theta)\big)\ge\mathcal L(q^{(t)},\theta), $$ বিশেষত $\theta=\theta^{(t+1)}$-এ: $\log p(x\mid\theta^{(t+1)})\ge\mathcal L(q^{(t)},\theta^{(t+1)})$। এবার $(\star\star)$ ও $(\star)$ জুড়ে: $$ \log p(x\mid\theta^{(t+1)})\;\ge\;\mathcal L(q^{(t)},\theta^{(t+1)})\;\overset{(\star\star)}{\ge}\;\mathcal L(q^{(t)},\theta^{(t)})\;\overset{(\star)}{=}\;\log p(x\mid\theta^{(t)}). $$ তাই $$ \boxed{\log p(x\mid\theta^{(t+1)})\ge\log p(x\mid\theta^{(t)})}\quad\blacksquare $$ — EM প্রতি iteration-এ log-likelihood একঘেয়ে বাড়ায় বা স্থির রাখে (কখনো কমায় না)। মূল উপাদান দুটো: E-step bound tight করে ($\star$), M-step bound বাড়ায় ($\star\star$), আর $\mathrm{KL}\ge0$ (Jensen) সেতু গড়ে। [টীকা: এটা monotone উন্নতি নিশ্চিত করে, কিন্তু global maximum নয় — EM local optimum-এ আটকাতে পারে, তাই একাধিক random restart (n_init) দরকার।]

ঘ · কোডিং (Python)¶

সমাধান ১৩ (★★★) — পূর্ণ pipeline¶

# 06-07 — Density Estimation, EM & GMM: GMM fit, responsibility, ARI vs k-means, BIC sweep
# seed-সূচক 20260619 — 3-component 2D Gaussian mixture, n=600, weights [0.40,0.35,0.25], tilted covs
import numpy as np
from sklearn.mixture import GaussianMixture
from sklearn.cluster import KMeans
from sklearn.metrics import adjusted_rand_score

SEED = 20260619
rng = np.random.RandomState(SEED)

# ----- canonical dataset: 3-component 2D Gaussian mixture, unequal/tilted covariances -----
n = 600
weights = [0.40, 0.35, 0.25]
means = [np.array([0.0, 0.0]),
         np.array([4.0, 4.0]),
         np.array([-3.0, 3.0])]
covs = [np.array([[1.0,  0.6], [ 0.6, 1.0]]),   # tilted (+ correlation)
        np.array([[1.5, -0.7], [-0.7, 0.8]]),   # tilted (- correlation), unequal scale
        np.array([[0.7,  0.2], [ 0.2, 2.0]])]   # elongated along y
counts = rng.multinomial(n, weights)
X_parts, y_parts = [], []
for k, c in enumerate(counts):
    X_parts.append(rng.multivariate_normal(means[k], covs[k], size=c))
    y_parts += [k] * c
X = np.vstack(X_parts)
y_true = np.array(y_parts)

# ----- 1. GMM fit (K=3) -----
gm = GaussianMixture(n_components=3, covariance_type='full',
                     random_state=SEED, n_init=10).fit(X)
print("=== GMM (K=3) ===")
print("BIC            :", round(gm.bic(X), 1))          # canonical 4828.8
print("per-sample LL  :", round(gm.score(X), 3))        # canonical -3.933
print("weights (↓)    :", sorted(np.round(gm.weights_, 3), reverse=True))  # ~[0.405,0.349,0.246]

# ----- 2. responsibilities + most-ambiguous point -----
resp = gm.predict_proba(X)                 # shape (600, 3) — soft assignment
amb = np.argmin(resp.max(axis=1))          # বিন্দু যার সর্বোচ্চ responsibility সবচেয়ে কম
print("\n=== responsibilities ===")
print("most ambiguous point γ:", np.round(resp[amb], 3))  # canonical ~[0.864,0.002,0.133]

# ----- 3. ARI: GMM (soft→hard) vs k-means (hard) -----
gm_labels = gm.predict(X)
km = KMeans(n_clusters=3, random_state=SEED, n_init=10).fit(X)
ari_gmm = adjusted_rand_score(y_true, gm_labels)
ari_km  = adjusted_rand_score(y_true, km.labels_)
print("\n=== ARI vs ground truth ===")
print("GMM  ARI:", round(ari_gmm, 3))      # canonical 0.97
print("kmeans ARI:", round(ari_km, 3))     # canonical 0.914
print("=> GMM (soft/elliptical) > k-means (hard/spherical)" if ari_gmm > ari_km else "")

# ----- 4. BIC sweep over K -----
print("\n=== BIC by K (lower = better) ===")
bics = {}
for K in range(1, 7):
    g = GaussianMixture(n_components=K, covariance_type='full',
                        random_state=SEED, n_init=10).fit(X)
    bics[K] = g.bic(X)
    print(f"K={K}: BIC = {round(bics[K], 1)}")
best_K = min(bics, key=bics.get)
print("best K (min BIC):", best_K)         # canonical K=3

প্রত্যাশিত আউটপুট (canonical):

=== GMM (K=3) ===
BIC            : 4828.8
per-sample LL  : -3.933
weights (↓)    : [0.405, 0.349, 0.246]

=== responsibilities ===
most ambiguous point γ: [0.864 0.002 0.133]

=== ARI vs ground truth ===
GMM  ARI: 0.97
kmeans ARI: 0.914
=> GMM (soft/elliptical) > k-means (hard/spherical)

=== BIC by K (lower = better) ===
K=1: BIC = 5626.7
K=2: BIC = 5111.7
K=3: BIC = 4828.8
K=4: BIC = 4857.5
K=5: BIC = 4890.9
K=6: BIC = 4925.9
best K (min BIC): 3

পাঠ। (১) GMM($3$) BIC $4828.8$, per-sample log-lik $-3.933$, আনুমানিক weights $[0.405,0.349,0.246]$ — প্রকৃত $[0.40,0.35,0.25]$-এর কাছাকাছি, fit ভালো। (২) সবচেয়ে অনিশ্চিত বিন্দুর responsibility $[0.864,0.002,0.133]$ — soft assignment-এর জীবন্ত দৃষ্টান্ত (দুই component-এর সীমানা)। (৩) GMM-ARI $0.97$ $>$ k-means-ARI $0.914$ — tilted/elliptical cluster-এ soft/elliptical GMM hard/spherical k-means-কে হারায়। (৪) BIC $K=3$-এ সর্বনিম্ন (V-আকার: $1\to3$ পড়ে, $3\to6$ ওঠে), ঠিক প্রকৃত cluster-সংখ্যা। (দ্রষ্টব্য: n_init=10 — EM non-convex, একাধিক random start না দিলে local optimum-এ আটকে আলাদা সংখ্যা আসতে পারে; data-recipe/সংস্করণ-ভেদে শেষ অঙ্কে সামান্য হেরফের সম্ভব, কিন্তু গল্প অপরিবর্তিত।)

সমাধান ১৪ (★★★) — নিজের হাতে EM (scratch, numpy)¶

# 06-07 — scratch 2-component 1D Gaussian-mixture EM (sklearn ছাড়া, শুধু numpy)
import numpy as np

x = np.array([1.0, 1.5, 5.0, 6.0, 6.5])   # দুটো স্পষ্ট গুচ্ছ: {1.0,1.5} ও {5.0,6.0,6.5}
n = len(x)

# প্রাথমিক parameter
pi  = np.array([0.5, 0.5])
mu  = np.array([0.0, 6.0])
var = np.array([1.0, 1.0])

def gauss(x, m, v):                        # 1D Gaussian density N(x; m, v)  (v = variance)
    return np.exp(-(x - m)**2 / (2*v)) / np.sqrt(2*np.pi*v)

def log_likelihood(x, pi, mu, var):
    p = sum(pi[k]*gauss(x, mu[k], var[k]) for k in range(2))
    return np.sum(np.log(p))

print(f"iter 0 : loglik = {log_likelihood(x, pi, mu, var):.4f}")
for it in range(1, 6):
    # ---- E-step: responsibilities γ (shape n×2) ----
    g = np.zeros((n, 2))
    for k in range(2):
        g[:, k] = pi[k] * gauss(x, mu[k], var[k])
    g /= g.sum(axis=1, keepdims=True)      # row-wise normalize ⇒ Σ_k γ_ik = 1

    # ---- M-step: weighted updates ----
    Nk  = g.sum(axis=0)                     # effective counts
    pi  = Nk / n
    mu  = (g * x[:, None]).sum(axis=0) / Nk
    var = (g * (x[:, None] - mu)**2).sum(axis=0) / Nk

    ll = log_likelihood(x, pi, mu, var)
    print(f"iter {it} : loglik = {ll:.4f}  mu = {np.round(mu,3)}  pi = {np.round(pi,3)}")

প্রত্যাশিত আউটপুট:

iter 0 : loglik = -10.3103
iter 1 : loglik = -6.2705  mu = [1.25  5.833]  pi = [0.4 0.6]
iter 2 : loglik = -6.2705  mu = [1.25  5.833]  pi = [0.4 0.6]
iter 3 : loglik = -6.2705  mu = [1.25  5.833]  pi = [0.4 0.6]
iter 4 : loglik = -6.2705  mu = [1.25  5.833]  pi = [0.4 0.6]
iter 5 : loglik = -6.2705  mu = [1.25  5.833]  pi = [0.4 0.6]

পাঠ। (ক) E-step প্রতিটি বিন্দুর responsibility $\gamma_{ik}$ গণনা করে ($\sum_k\gamma_{ik}=1$, row-normalize)। (খ) M-step $N_k=\sum_i\gamma_{ik}$ (effective count), $\pi_k=N_k/n$, weighted-mean $\mu_k$, weighted-variance $\sigma_k^2$ আপডেট করে — §গ-এর derive-করা সূত্র (সমাধান ১১)। (গ) log-likelihood একঘেয়ে বাড়ছে ($-10.31\to-6.27$, তারপর স্থির — এই পরিষ্কার data এক পূর্ণ sweep-এই অভিসৃত), ঠিক সমাধান ১২-এর তত্ত্ব অনুযায়ী ($\ell^{(t+1)}\ge\ell^{(t)}$, কখনো কমে না)। দুটো cluster সঠিকভাবে আলাদা হয়: $\mu\to(1.25,\,5.833)$ — গুচ্ছ $\{1.0,1.5\}$-এর গড় $1.25$ আর $\{5.0,6.0,6.5\}$-এর গড় $5.833$, weights $[0.4,0.6]$ ($2$ বনাম $3$ বিন্দু)। এই খেলনা-EM থেকেই §১৩-এর sklearn GMM-এর অভ্যন্তরীণ যন্ত্রটা পরিষ্ক

দিক	KDE (kernel density)	mixture model (GMM)
(ক) component-সংখ্যা	\(n\) — প্রতি data-বিন্দুতে এক kernel; data বাড়লে সংখ্যাও বাড়ে ⇒ nonparametric	ধ্রুবক \(K\) (অল্প) — data-আকার-নিরপেক্ষ ⇒ parametric
(খ) প্রতি মূল্যায়নে term	\(n\)টি (পুরো data সঞ্চয়/স্ক্যান করতে হয়)	মাত্র \(K\)টি Gaussian — বড় \(n\)-এ অনেক সাশ্রয়ী ও compact