Mai jos regăsiți subiectele pentru proiectul de laborator pentru cursul de Statistică.
Fiecare echipa va trimite prin liderul echipei un singur e-mail la adresa simona.cojocea@fmi.unibuc.ro ce va conține o arhivă cu 2 componente: codul
R(comentat!) și documentația proiectului, în format .docx sau .pdf.Documentația proiectului trebuie să conțină, pe prima pagină, numele membrilor echipei, liderul echipei și grupa din care face parte fiecare. Pentru fiecare exercițiu în parte documentația trebuie să conțină:
- Calculele matematice solicitate (pot fi tehnoredactate sau scrise de mână și scanate, însă dacă apelați la ultima variantă vă rog să vă asigurați că ați scris citeț, fără ștersături sau zone scanate deficitar)
- Codul
Rcomentat - Graficele realizate (acolo unde au fost cerute)
- Comentariile și concluziile voastre
Pentru obținerea punctajului maxim trebuie să rezolvați toate subiectele. Dacă odată cu rezolvarea lor, identificați și rezolvați unele cerințe suplimentare care sunt relevante pentru subiectul respectiv puteți obține un bonus de până la 10p, fără ca nota finală de laborator să poată depăși 50p.
Data de predare a proiectului este 4 februarie 2024 ora 22:00.
Subiecte
Exercițiul 1 Construiți funcții în R care să implementeze algoritmi de simulare pentru următoarele situații descrise:
\(X \sim \operatorname{Cauchy}\left(x_0, \gamma\right)\)
\(X \sim\left(\begin{array}{ccc}1 & 2 & 3 \\ \frac{1}{2} & \frac{1}{3} & \frac{1}{6}\end{array}\right)\)
\(X \sim\left(\begin{array}{cc}x_1 & x_2 \\ 1-p & p\end{array}\right)\) unde \(x_1\), \(x_2\) și \(p\) sunt transmiși funcției ca parametri
extragerea unui număr întreg cuprins între \(a\) și \(b\) într-o manieră perturbată aleator față de repartiția uniformă (de exemplu, se extrage ora curentă a sistemului, se face o prelucrare asupra acesteia pentru a obține un număr natural \(x\) care va reprezenta poziția din șirul numerelor naturale cuprinse între \(a\) și \(b\) căreia i se va asocia o probabilitate mai mare de apariție)
În continuare, veți aplica acești algoritmi de simulare pentru a rezolva următoarea problemă:
Despre linia 501 STB se cunosc:
în medie călătoresc \(y\) călători pe zi (medie raportată la datele dintr-o lună pentru un singur tramvai), unde \(y\) se obține apelând funcția de la d)
numărul minim înregistrat în perioada studiată este de \(x_{min}\) iar numărul maxim este \(x_{max}\) călători pe zi (valori ce se obțin folosind funcția de la d))
zilele dintr-o lună sunt clasificate ca fiind lejere (mai puțin de 350 de călători pe zi), normale (între 351 și 670 călători pe zi), aglomerate (peste 671 călători pe zi). Folosiți funcția de la b) pentru a simula câte zile lejere, normale și respectiv aglomerate are fiecare lună în parte.
prețul unui bilet este 3 lei și, în medie, \(x\%\) din pasagerii care nu au abonament plătesc biletul la utilizarea tramvaiului, unde \(x\) este un număr între 1 și 99 generat astfel: cu algoritmul de la a) se generează un număr v căruia i se extrage partea fracționară și i se atribuie parametrului p de la c), iar \(x_1\) și respectiv \(x_2\) sunt extrase uniform din mulțimea numerelor naturale de la 1 la 99. Valoarea medie va fi rezultatul simulării de la c) în condițiile descrise.
prețul unui abonament este 70 lei pe lună și, în medie, \(x\%\) din pasagerii care călătoresc cu tramvaiul îl achiziționează, unde \(x\) este un număr între 1 și 99 generat cu algoritmul de la d) folosind ca element perturbator un număr z obținut ca modulul părții întregi a unui număr obținut prin apelul funcției de la a) cu parametrii 5 și 2.
Se cere:
Generați, prin simulare, valori care să reprezinte numărul de călători dintr-o zi, pentru fiecare zi a lunii decembrie 2024, respectând restricțiile de mai sus și stocați valorile obținute într-un vector. Construiți histograma acestor valori.
Repetați procedeul de la a) pentru fiecare luna a anului 2024 și centralizați rezultatele empirice într-un
data.framecare să conțină, pentru fiecare lună valorile medii, minime și maxime de călători, precum și procentul de zile lejere, normale și respectiv aglomerate înregistrate.Completați
data.frame-ul de la b) cu simularea numărului de pasageri cu abonament, numărului de pasageri care plătesc bilet și respectiv numărului de pasageri care nu plătesc bilet. Determinați pentru fiecare lună în parte veniturile provenite din bilete și abonamente și respectiv, veniturile nerealizate prin neplata biletului de unii dintre pasageri. Organizați informația într-o manieră ușor de vizualizat.Un tramvai de pe linia lui 501 face 14 trasee complete în timpul programului de lucru dintr-o zi. De două ori pe zi un controlor se urcă în unul din tramvaie și solicită prezentarea biletelor de călătorie unui număr de pasageri, aleși în mod aleator, după următorul algoritm:
- dacă e o zi lejeră, verifică în mod aleator între 2 și 11 persoane, dar se oprește din verificare dacă a amendat deja 3 persoane.
- dacă e zi normală, verifică în mod aleator un număr de persoane până reusește să amendeze 5 persoane (sau a verificat pe toata lumea prezentă în tramvai între 2 stații)
- dacă e o zi aglomerată, verifică în mod aleator între 3 și 5 persoane și se oprește din verificare după prima amendă
Știind că amenda este 50 lei, determinați, în urma simulării, pentru fiecare zi a fiecărei luni din anul 2025, câți bani se strâng din aplicarea unor amenzi. Comparați această sumă cu pierderea realizată prin neplata biletelor și stabiliți în câte zile dintr-o lună (în medie) sumele obținute din amenzi depășesc pierderea prin neplata biletelor, considerând și faptul că pentru fiecare din cele 2 verificări zilnice există un cost asociat controlorului de 214 lei.
- Studiați prin simulare, oportunitatea de a introduce un al treilea control pe zi (în condițiile menționate anterior, la care adaugăm informația ca, în medie, 30% din amenzile colectate la fiecare din controale nu sunt raportate oficial ci sunt păstrate de controlor).
Pentru exercițiul următor vom introduce noțiunea de familie exponențială:
Definiția 1 (Familia exponențială) O variabilă aleatoare \(X\) face parte din familia exponențială s-dimensională dacă densitatea/funcția de masă poate fi scrisă sub forma:
\[ p(x ; \theta)=h(x) \cdot \exp \left(\sum_{i=1}^s \eta_i(\theta) \cdot T_i(x)-A(\theta)\right) \]
unde \(\eta_i\) și \(A\) sunt funcții reale de \(\theta=\left(\theta_1, \theta_2, \ldots \theta_s\right)\), \(T_i\) reprezintă statistici suficiente, iar \(h\) este o funcție pozitivă de \(x\).
Funcția \(A(\theta)\) se numește constantă de log-normalizare (rolul ei este acela de a face ca să îndeplinească acele condiții necesare pentru a fi o funcție de densitate de probabilitate/funcție de masă, după caz).
Exemplu: Fie \(X\sim\mathcal{N}(0,1)\). Vom arăta că aceasta face parte din familia exponențială 2-dimensională. Avem
\[ \begin{aligned} & f\left(x ; \mu, \sigma^2\right)=\frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left(\frac{(x-\mu)^2}{2 \sigma^2}\right)\\ &=\frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left(\frac{\mu}{\sigma^2} x \frac{1}{2 \sigma^2} x^2 \frac{\mu^2}{2 \sigma^2}\right) \\ & =\frac{1}{\sqrt{2 \pi}} \exp \left(\frac{\mu}{\sigma^2} x-\frac{1}{2 \sigma^2} x^2-\frac{\mu^2}{2 \sigma^2}-\ln (\sigma)\right) \end{aligned} \]
de unde identificăm următoarele relații:
\[ \begin{aligned} & h(x)=\frac{1}{\sqrt{2 \pi}} \\ & \eta_1\left(\theta_1, \theta_2\right)=\frac{\mu}{\sigma^2},\quad \eta_2\left(\theta_1, \theta_2\right)=-\frac{1}{2 \sigma^2},\quad \theta_1=\mu,\quad \theta_2=\sigma^2 \\ & A(\theta)=\frac{\mu^2}{2 \sigma^2}+\ln (\sigma) \\ & T_1(x)=x, T_2(x)=x^2 \end{aligned} \]
Așadar, repartiția normală (cu parametrii media \(\mu\) și dispersia \(\sigma^2\)) face parte din familia exponențială 2-dimensională.
Exercițiul 2 Ne plasăm în contextul definiției de mai sus.
Verificați dacă următoarele repartiții fac parte din familia exponențială:
- \(\mathcal{B}(3,p)\)
- \(\mathcal{B}(n,p)\)
- \(\mathrm{Geom}(p)\)
- \(\mathrm{Pois}(\lambda)\)
- \(\Gamma(\alpha,\beta)\)
- \(\mathrm{Beta}(\alpha, \beta)\)
- \(\chi^2(\nu)\)
Ilustrați grafic în
Rdensitățile/funcțiile de masă (după caz) ale repartițiilor de mai sus pentru 4 parametri particulari, la alegere, în cadrul aceluiași sistem de axe ortogonale (fiecare repartiție va avea însă reprezentări grafice distincte de celelalte repartiții).Construiți funcția de log-verosimilitate (logL) pentru familia exponențială și demonstrați că aceasta este concavă (construiți matricea hessiană și arătați că este negativ definită).
Particularizați forma funcției de log-verosimilitate de la 2) pentru repartițiile de la 1) ce fac parte din clasa exponențială cu un parametru. Reprezentați grafic (în
R) aceste funcții și găsiți punctul lor de maxim (folosiți funcțiaoptimize).Particularizați forma funcției de log-verosimilitate de la 2) pentru repartițiile de la 1) ce fac parte din clasa exponențială cu doi parametri. Fixați, pe rând, unul din parametri (alegeți o valoare particulară după cum doriți), reprezentați grafic (în
R) aceste funcții în raport cu celălalt parametru și găsiți punctul lor de maxim (folosiți funcțiaoptimize).Calculați MIRC pentru familia exponențială cu un parametru și particularizați valoarea acesteia pentru repartițiile de la 1) ce fac parte din familia exponențială cu un parametru.
Construiți în
Ro funcție care afișează MIRC pentru o repartiție selectată din 8 disponibile (alegeți voi aceste repartiții).Pentru familia exponențială cu un parametru și respectiv cu doi parametri calculați estimatorul de verosimilitate maximă și respectiv estimatorul dat de metoda momentelor. Ce legătura există între aceștia?
Ilustrați în
R, pentru un eșantion de volum \(1000\) generat de voi în prealabil - pentru toate repartițiile de la 1) ce fac parte din familia exponențială - faptul că estimațiile obținute în baza celor 2 metode de estimare (metoda verosimilității maxime și metoda momentelor) pentru eșantionul respectiv sunt foarte apropiate de valoarea adevărata a parametrului de interes.
Exercițiul 3 Inegalitatea Berry-Esseen este un rezultat celebru cu ajutorul căruia putem determina acuratețea aproximării pe care o realizează TLC. Cerând în plus față de condițiile din TLC ca \(\mathbb{E}[|X_i|^3]<\infty\), avem că:
\[ \sup_{x}\left|\mathbb{P}(Z_n\leq x) - \Phi(x)\right|\leq \frac{33}{4}\frac{\mathbb{E}[|X_1-\mu|^3]}{\sigma^3\sqrt{n}} \]
unde \(X_1,\ldots,X_n\) i.i.d., \(\Phi(x)\) este funcția de repartiție a normalei standard, \(\mu=\mathbb{E}[X_1]\), \(\sigma = \sqrt{Var{X_1}}\) și \(Z_n=\frac{\sqrt{n}\left(\bar{X}_n-\mu\right)}{\sigma}\).
Calculați \(\mathbb{P}(Z_n\leq x)\) pentru repartițiile: Binomială, Geometrică, Poisson, Uniformă pe caz discret și respectiv continuu, Exponențială, Gamma și respectiv Beta.
Reprezentați grafic funcțiile obținute la punctul 1.
Folosind funcția
optimizeaproximați pentru repartițiile de mai sus: \(\sup_{x}\left|\mathbb{P}(Z_n\leq x) - \Phi(x)\right|\).Construiți câte o funcție în
Rcare să calculeze \(\mathbb{E}[X]\) și respectiv \(Var(X)\), unde tipul repartiției v.a. \(X\) este transmis fie printr-o denumire, fie prin funcția de masă în cazul discret, fie prin funcția densitate de probabilitate în cazul continuu.Construiți o funcție în
Rcare să calculeze \(\mathbb{E}[|X_1-\mu|^3]\)Construiți un
data.frameînRcare să conțină marginile date de inegalitatea Berry-Esseen pentru repartițiile: Binomială, Geometrică, Poisson, Uniformă pe caz discret și respectiv Exponențială, Gamma și Beta în cazul continuu, pentru un volum al eșantionului \(n\in\{30, 100, 1000\}\). Alegerea valorilor parametrilor repartițiilor ilustrate vă revine vouă.Ilustrați grafic pentru \(n\in\{30, 100, 1000\}\) pentru un interval pe care îl considerați relevant evoluția diferenței \(\mathbb{P}(Z_n\leq x) - \Phi(x)\) pentru fiecare din repartițiile de mai sus.
Construiți o funcție în
Rcare să calculeze marginea dată de inegalitatea Berry-Esseen pentru o v.a. pentru care se cunoaște funcția de masă/funcția densitate de probabilitate.