Appendice Y — Probabilità

# Standard library imports
import os

# Third-party imports
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import arviz as az
import scipy.stats as stats
from scipy.special import expit  # Funzione logistica
import math
from cmdstanpy import cmdstan_path, CmdStanModel

# Configuration
seed = sum(map(ord, "stan_poisson_regression"))
rng = np.random.default_rng(seed=seed)
az.style.use("arviz-darkgrid")
%config InlineBackend.figure_format = "retina"

# Define directories
home_directory = os.path.expanduser("~")
project_directory = f"{home_directory}/_repositories/psicometria"

# Print project directory to verify
print(f"Project directory: {project_directory}")

Project directory: /Users/corradocaudek/_repositories/psicometria

Capitolo 27

Esercizio 27.1

Per calcolare questa probabilità in maniera analitica, utilizziamo la seguente uguaglianza:

$P (almeno 2 psicologi clinici) = 1 - P (nessun psicologo clinico) - P (1 psicologo clinico) .$

Il numero totale di modi per selezionare 5 persone dal gruppo di 20 è dato da:

$(\binom{20}{5}) = \frac{20!}{5! (15!)} = 15, 504.$

Il numero di modi per avere nessun psicologo clinico nella commissione (ovvero, selezionare solo psicologi del lavoro) è:

$(\binom{10}{0}) \times (\binom{10}{5}) = 1 \times 252 = 252.$

Quindi, la probabilità di avere nessun psicologo clinico è:

$P (nessun psicologo clinico) = \frac{252}{15, 504} \approx 0.016 .$

Il numero di modi per avere esattamente 1 psicologo clinico nella commissione è:

$(\binom{10}{1}) \times (\binom{10}{4}) = 10 \times 210 = 2, 100.$

Quindi, la probabilità di avere esattamente 1 psicologo clinico è:

$P (1 psicologo clinico) = \frac{2, 100}{15, 504} \approx 0.135 .$

La probabilità di avere almeno 2 psicologi clinici nella commissione è quindi:

$\begin{aligned} P (almeno 2 psicologi clinici) & = 1 - P (nessun psicologo clinico) - P (1 psicologo clinico) \\ = 1 - 0.016 - 0.135 \\ = 0.848 . \end{aligned}$

Quindi, la probabilità che almeno 2 psicologi clinici siano nella commissione è circa 0.848.

# Funzione per calcolare le combinazioni
def nCk(n, k):
    return math.factorial(n) // (math.factorial(k) * math.factorial(n - k))


# Calcolo delle probabilità per il problema della commissione
total_ways = nCk(20, 5)
no_clinical = nCk(10, 0) * nCk(10, 5)
one_clinical = nCk(10, 1) * nCk(10, 4)

p_no_clinical = no_clinical / total_ways
p_one_clinical = one_clinical / total_ways

p_at_least_two_clinical = 1 - p_no_clinical - p_one_clinical

print(f"Probabilità di almeno 2 psicologi clinici: {p_at_least_two_clinical:.3f}")

Probabilità di almeno 2 psicologi clinici: 0.848

In maniera più intuitiva, possiamo risolvere il problema con una simulazione Monte Carlo.

import random

# Numero di simulazioni
simulations = 1000000

# Numero di successi (almeno 2 psicologi clinici nella commissione)
success_count = 0

# Creiamo una lista che rappresenta il gruppo di 20 persone
# 1 rappresenta un psicologo clinico, 0 rappresenta un psicologo del lavoro
group = [1] * 10 + [0] * 10

# Simulazione Monte Carlo
for _ in range(simulations):
    # Estrai casualmente 5 persone dal gruppo
    committee = random.sample(group, 5)

    # Conta quanti psicologi clinici ci sono nella commissione
    num_clinical_psychologists = sum(committee)

    # Verifica se ci sono almeno 2 psicologi clinici
    if num_clinical_psychologists >= 2:
        success_count += 1

# Calcola la probabilità
probability = success_count / simulations

# Mostra il risultato
print(
    f"La probabilità che almeno 2 psicologi clinici siano nella commissione è: {probability:.4f}"
)

La probabilità che almeno 2 psicologi clinici siano nella commissione è: 0.8482

Capitolo 39

Esercizio 39.1

Per calcolare le deviazioni standard delle distribuzioni gaussiane date le percentuali di studenti che ottengono meno di 18, possiamo utilizzare le proprietà della distribuzione normale e i quantili della distribuzione normale standard (distribuzione normale con media 0 e deviazione standard 1).

Le distribuzioni normali hanno la proprietà che possiamo trasformare qualsiasi valore $X$ della distribuzione $N (μ, σ)$ nella distribuzione normale standard $N (0, 1)$ tramite la formula:

$Z = \frac{X - μ}{σ},$

dove $Z$ è il quantile standardizzato.

Per trovare il valore di $σ$ dato un certo percentile, utilizziamo l’inverso della funzione di distribuzione cumulativa (CDF) della distribuzione normale standard. Per un dato percentile $p$ , $z_{p}$ è tale che:

$p = P (Z \leq z_{p})$

Quindi possiamo trovare $σ$ risolvendo per $σ$ nella formula:

$z_{p} = \frac{X - μ}{σ},$

$σ = \frac{X - μ}{z_{p}},$

dove:

$X$ è il punteggio di soglia (18 in questo caso).
$μ$ è la media della distribuzione.
$z_{p}$ è il quantile della distribuzione normale standard per il percentile $p$ .

I quantili della distribuzione normale standard per i percentili desiderati sono:

Per il 15%, il quantile è $z_{0.15} \approx - 1.036$ .
Per il 10%, il quantile è $z_{0.10} \approx - 1.281$ .
Per il 5%, il quantile è $z_{0.05} \approx - 1.645$ .

Prima Prova

Media: $μ = 24$
Percentuale che ottiene meno di 18: 15%
Quantile: $z_{0.15} = - 1.036$
Soglia: $X = 18$

$σ_{1} = \frac{24 - 18}{1.036} \approx 5.79$

Seconda Prova

Media: $μ = 25$
Percentuale che ottiene meno di 18: 10%
Quantile: $z_{0.10} = - 1.281$
Soglia: $X = 18$

$σ_{2} = \frac{25 - 18}{1.281} \approx 5.46$

Terza Prova

Media: $μ = 26$
Percentuale che ottiene meno di 18: 5%
Quantile: $z_{0.05} = - 1.645$
Soglia: $X = 18$

$σ_{3} = \frac{26 - 18}{1.645} \approx 4.86$

# Funzione per calcolare la deviazione standard data la media, la soglia e il quantile
def calculate_std(mean, threshold, quantile):
    return abs((mean - threshold) / quantile)


# Parametri delle distribuzioni gaussiane per le tre prove
mean_test1 = 24
std_test1 = calculate_std(mean_test1, 18, -1.036)
mean_test2 = 25
std_test2 = calculate_std(mean_test2, 18, -1.281)
mean_test3 = 26
std_test3 = calculate_std(mean_test3, 18, -1.645)

# Numero di studenti
n_students = 220

# Percentuale di studenti che non fa le prove
drop_test1 = 0.10
drop_test2 = 0.05

# Seed per il generatore di numeri casuali basato sulla stringa "simulation"
seed = sum(map(ord, "simulation"))
rng = np.random.default_rng(seed=seed)

# Generazione dei voti per le tre prove
# Genera i voti solo per gli studenti che partecipano alla prova
test1_scores = np.where(
    rng.random(n_students) > drop_test1,
    rng.normal(mean_test1, std_test1, n_students),
    np.nan,
)
test2_scores = np.where(
    rng.random(n_students) > drop_test2,
    rng.normal(mean_test2, std_test2, n_students),
    np.nan,
)
test3_scores = rng.normal(mean_test3, std_test3, n_students)

# Calcola il voto finale solo per gli studenti che hanno partecipato a tutte e tre le prove
final_scores = np.nanmean(
    np.column_stack((test1_scores, test2_scores, test3_scores)), axis=1
)

# Filtra gli studenti che non hanno partecipato a tutte e tre le prove
valid_final_scores = final_scores[~np.isnan(final_scores)]

# Visualizzazione della distribuzione finale dei voti
plt.hist(valid_final_scores, bins=30, edgecolor="black")
plt.title("Distribuzione dei voti finali")
plt.xlabel("Voto finale")
plt.ylabel("Frequenza")
plt.show()

# Statistiche descrittive dei voti finali
mean_final_score = np.mean(valid_final_scores)
median_final_score = np.median(valid_final_scores)
std_final_score = np.std(valid_final_scores)

print(f"Media dei voti finali: {mean_final_score:.2f}")
print(f"Mediana dei voti finali: {median_final_score:.2f}")
print(f"Deviazione standard dei voti finali: {std_final_score:.2f}")