import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats
import math


VARIANT = 2
VARIABLE_NAME = "Вага, кг"

UNIFORM_MIN = 40
UNIFORM_MAX = 90

N_OBS = 100

RANDOM_SEED = 42

def generate_uniform_data(min_val: float, max_val: float, n: int, seed: int = None) -> np.ndarray:
    rng = np.random.default_rng(seed)
    raw = rng.uniform(min_val, max_val, n)
    rounded = np.round(raw, 0).astype(int)
    return rounded


def stats_by_formulas(x: np.ndarray) -> dict:
    x = np.asarray(x, dtype=float)
    n = x.size
    x_sorted = np.sort(x)
    x_min = x_sorted[0]
    x_max = x_sorted[-1]
    R = x_max - x_min

    x_mean = x.sum() / n

    if n % 2 == 1:
        median = x_sorted[n // 2]
    else:
        median = 0.5 * (x_sorted[n // 2 - 1] + x_sorted[n // 2])

    values, counts = np.unique(x, return_counts=True)
    mode = values[np.argmax(counts)]

    d = np.mean(np.abs(x - x_mean))

    var_pop = np.mean((x - x_mean) ** 2)
    sigma_pop = np.sqrt(var_pop)

    var_unbiased = var_pop * n / (n - 1)
    sigma_unbiased = np.sqrt(var_unbiased)

    V_R = (R / x_mean) * 100.0
    V_d = (d / x_mean) * 100.0
    V_sigma = (sigma_pop / x_mean) * 100.0

    centered = x - x_mean
    A = (centered ** 3).sum() / (n * (sigma_pop ** 3)) if sigma_pop != 0 else np.nan

    E = (centered ** 4).sum() / (n * (sigma_pop ** 4)) - 3 if sigma_pop != 0 else np.nan

    return {
        "n": n,
        "min": x_min,
        "max": x_max,
        "R": R,
        "mean": x_mean,
        "median": median,
        "mode": mode,
        "d": d,
        "var_pop": var_pop,
        "sigma_pop": sigma_pop,
        "var_unbiased": var_unbiased,
        "sigma_unbiased": sigma_unbiased,
        "V_R": V_R,
        "V_d": V_d,
        "V_sigma": V_sigma,
        "A": A,
        "E": E,
    }


def stats_by_functions(x: np.ndarray) -> dict:
    x = np.asarray(x, dtype=float)
    n = x.size
    x_min = np.min(x)
    x_max = np.max(x)
    R = x_max - x_min

    x_mean = x.mean()
    median = np.median(x)
    mode = float(pd.Series(x).mode().iloc[0])

    d = np.mean(np.abs(x - x_mean))

    var_pop = np.var(x, ddof=0)
    sigma_pop = np.std(x, ddof=0)
    var_unbiased = np.var(x, ddof=1)
    sigma_unbiased = np.std(x, ddof=1)

    V_R = (R / x_mean) * 100.0
    V_d = (d / x_mean) * 100.0
    V_sigma = (sigma_pop / x_mean) * 100.0

    A = stats.skew(x, bias=False)
    E = stats.kurtosis(x, fisher=True, bias=False)

    return {
        "n": n,
        "min": x_min,
        "max": x_max,
        "R": R,
        "mean": x_mean,
        "median": median,
        "mode": mode,
        "d": d,
        "var_pop": var_pop,
        "sigma_pop": sigma_pop,
        "var_unbiased": var_unbiased,
        "sigma_unbiased": sigma_unbiased,
        "V_R": V_R,
        "V_d": V_d,
        "V_sigma": V_sigma,
        "A": A,
        "E": E,
    }


def stats_by_addin(x: np.ndarray) -> dict:
    x = np.asarray(x, dtype=float)
    desc = stats.describe(x)
    n = desc.nobs
    x_min, x_max = desc.minmax
    R = x_max - x_min
    x_mean = desc.mean

    var_unbiased = float(desc.variance)
    sigma_unbiased = float(np.sqrt(var_unbiased))

    var_pop = var_unbiased * (n - 1) / n
    sigma_pop = float(np.sqrt(var_pop))

    A = float(desc.skewness)
    E = float(desc.kurtosis)

    median = float(np.median(x))
    mode = float(pd.Series(x).mode().iloc[0])
    d = float(np.mean(np.abs(x - x_mean)))

    V_R = (R / x_mean) * 100.0
    V_d = (d / x_mean) * 100.0
    V_sigma = (sigma_pop / x_mean) * 100.0

    return {
        "n": n,
        "min": x_min,
        "max": x_max,
        "R": R,
        "mean": x_mean,
        "median": median,
        "mode": mode,
        "d": d,
        "var_pop": var_pop,
        "sigma_pop": sigma_pop,
        "var_unbiased": var_unbiased,
        "sigma_unbiased": sigma_unbiased,
        "V_R": V_R,
        "V_d": V_d,
        "V_sigma": V_sigma,
        "A": A,
        "E": E,
    }


from typing import Tuple

def check_homogeneity(V_sigma: float) -> str:
    if V_sigma <= 30.0:
        return f"Вибірка однорідна (V_sigma = {V_sigma:.2f}% ≤ 30%)"
    else:
        return f"Вибірка НЕ є однорідною (V_sigma = {V_sigma:.2f}% > 30%)"

def critical_values_A_E(n: int) -> Tuple[float, float]:
    D_A = 6 * (n - 1) / ((n + 1) * (n + 3))
    D_E = 24 * n * (n - 2) * (n - 3) / (((n + 1) ** 2) * (n + 3) * (n + 5))
    return D_A, D_E

def check_normality(A: float, E: float, n: int) -> str:
    D_A, D_E = critical_values_A_E(n)
    sigma_A = math.sqrt(D_A)
    sigma_E = math.sqrt(D_E)
    crit_A = 3 * sigma_A
    crit_E = 5 * sigma_E

    cond_A = abs(A) <= crit_A
    cond_E = abs(E) <= crit_E

    lines = []
    lines.append(f"A = {A:.4f}, D(A) = {D_A:.4f}, 3·sqrt(D(A)) = {crit_A:.4f}")
    lines.append(f"E = {E:.4f}, D(E) = {D_E:.4f}, 5·sqrt(D(E)) = {crit_E:.4f}")
    if cond_A and cond_E:
        lines.append("Форма розподілу МОЖЕ вважатися близькою до нормальної (обидві нерівності виконуються).")
    else:
        lines.append("Форма розподілу НЕ є близькою до нормальної (хоча б одна з нерівностей порушена).")
    return "\n".join(lines)

def sturges_intervals(n: int) -> int:
    if n <= 1:
        return 1
    k = 1 + 3.322 * math.log10(n)
    return int(round(k))

def build_intervals(x: np.ndarray, k: int) -> np.ndarray:
    x = np.asarray(x, dtype=float)
    x_min = float(np.min(x))
    x_max = float(np.max(x))
    boundaries = np.linspace(x_min, x_max, k + 1)
    boundaries = np.round(boundaries, 1)
    return boundaries

def grouping_manual(x: np.ndarray, boundaries: np.ndarray) -> pd.DataFrame:
    x = np.asarray(x, dtype=float)
    n = x.size
    k = len(boundaries) - 1

    intervals_left = []
    intervals_right = []
    midpoints = []
    freqs = []

    for i in range(k):
        left = boundaries[i]
        right = boundaries[i + 1]

        if i < k - 1:
            mask = (x >= left) & (x < right)
        else:
            mask = (x >= left) & (x <= right)

        count = int(mask.sum())

        intervals_left.append(left)
        intervals_right.append(right)
        midpoints.append(0.5 * (left + right))
        freqs.append(count)

    freqs = np.array(freqs, dtype=int)
    rel_freq = freqs / n
    cum_freq_pct = np.cumsum(rel_freq) * 100.0

    df = pd.DataFrame({
        "L": intervals_left,
        "U": intervals_right,
        "midpoint": midpoints,
        "n_manual": freqs,
        "p_manual": rel_freq,
        "K_manual_pct": cum_freq_pct,
    })
    return df

def grouping_addin(x: np.ndarray, boundaries: np.ndarray, df_manual: pd.DataFrame) -> pd.DataFrame:
    x = np.asarray(x, dtype=float)
    k = len(boundaries) - 1

    freqs_addin, _ = np.histogram(x, bins=boundaries)

    n = x.size
    rel_freq_addin = freqs_addin / n
    cum_freq_addin_pct = np.cumsum(rel_freq_addin) * 100.0

    df = df_manual.copy()
    df["n_addin"] = freqs_addin
    df["p_addin"] = rel_freq_addin
    df["K_addin_pct"] = cum_freq_addin_pct
    return df

def theoretical_uniform_counts(boundaries: np.ndarray, n: int, a: float, b: float) -> np.ndarray:
    probs = (boundaries[1:] - boundaries[:-1]) / (b - a)
    counts = n * probs
    return counts

def plot_hist_addin(x: np.ndarray, boundaries: np.ndarray):
    plt.figure(figsize=(8, 5))
    plt.hist(x, bins=boundaries, edgecolor="black", alpha=0.7)
    plt.title("Гістограма (аналог статистичної надбудови)")
    plt.xlabel(VARIABLE_NAME)
    plt.ylabel("Частота")
    plt.grid(axis="y", alpha=0.3)

def plot_hist_polygon_cumulative(df_group_manual: pd.DataFrame):
    fig, ax1 = plt.subplots(figsize=(8, 5))

    L = df_group_manual["L"].values
    U = df_group_manual["U"].values
    mid = df_group_manual["midpoint"].values
    n_manual = df_group_manual["n_manual"].values
    K_manual_pct = df_group_manual["K_manual_pct"].values

    widths = U - L
    ax1.bar(mid, n_manual, width=widths * 0.9, edgecolor="black", alpha=0.6, label="Гістограма (ручне)")

    ax1.plot(mid, n_manual, "o-", label="Полігон")
    ax1.set_xlabel(VARIABLE_NAME)
    ax1.set_ylabel("Частота")
    ax1.set_title("Гістограма та полігон (ручне групування)")
    ax1.grid(alpha=0.3)

    ax2 = ax1.twinx()
    ax2.plot(mid, K_manual_pct, "s-", color="purple", label="Кумулята")
    ax2.set_ylabel("Накопичена частота, %")

    lines1, labels1 = ax1.get_legend_handles_labels()
    lines2, labels2 = ax2.get_legend_handles_labels()
    ax1.legend(lines1 + lines2, labels1 + labels2, loc="upper left")

def plot_empirical_vs_theoretical(df_group_manual: pd.DataFrame, theor_counts: np.ndarray):
    plt.figure(figsize=(8, 5))

    mid = df_group_manual["midpoint"].values
    n_manual = df_group_manual["n_manual"].values

    plt.bar(mid, n_manual, width=(df_group_manual["U"].values - df_group_manual["L"].values) * 0.9,
            edgecolor="black", alpha=0.6, label="Емпірична гістограма")

    plt.plot(mid, theor_counts, "o-", label="Теоретичний розподіл")

    plt.xlabel(VARIABLE_NAME)
    plt.ylabel("Частота")
    plt.title("Емпіричний vs теоретичний розподіл")
    plt.grid(alpha=0.3)
    plt.legend()

def run_lab(show_plots: bool = True):
    x = generate_uniform_data(UNIFORM_MIN, UNIFORM_MAX, N_OBS, seed=RANDOM_SEED)
    df = pd.DataFrame({"X": x})

    print(f"Практична робота 1-3. Варіант {VARIANT}. Випадкова величина: {VARIABLE_NAME}")
    print(f"Розподіл: рівномірний на [{UNIFORM_MIN}; {UNIFORM_MAX}], n = {N_OBS}")
    print("\nПерші 10 спостережень (округлені значення):")
    print(df.head(10))

    stats_formulas = stats_by_formulas(x)
    stats_funcs = stats_by_functions(x)
    stats_addin = stats_by_addin(x)

    all_keys = list(stats_formulas.keys())
    table = pd.DataFrame({
        "Показник": all_keys,
        "За формулами": [stats_formulas[k] for k in all_keys],
        "Стат. функції": [stats_funcs[k] for k in all_keys],
        "Надбудова": [stats_addin[k] for k in all_keys],
    })
    table.set_index("Показник", inplace=True)

    print("\n=== Показники варіації (три способи) ===")
    print(table)

    homog_text = check_homogeneity(stats_formulas["V_sigma"])
    print("\n=== Однорідність вибірки ===")
    print(homog_text)

    normality_text = check_normality(stats_formulas["A"], stats_formulas["E"], stats_formulas["n"])
    print("\n=== Перевірка близькості до нормального розподілу ===")
    print(normality_text)

    k = sturges_intervals(N_OBS)
    print(f"\nКількість інтервалів (формула Стерджеса): k = {k}")
    boundaries = build_intervals(x, k=k)
    df_group_manual = grouping_manual(x, boundaries)
    df_group = grouping_addin(x, boundaries, df_group_manual)

    print("\n=== Групування даних (ручне + надбудова) ===")
    print(df_group)

    theor_counts = theoretical_uniform_counts(boundaries, N_OBS, UNIFORM_MIN, UNIFORM_MAX)

    print("\nТеоретичні частоти для рівномірного розподілу:")
    print(theor_counts)

    if show_plots:
        plot_hist_addin(x, boundaries)
        plot_hist_polygon_cumulative(df_group_manual)
        plot_empirical_vs_theoretical(df_group_manual, theor_counts)
        plt.show()


run_lab(show_plots=True)





print(f"--- РОЗВ'ЯЗАННЯ КОНТРОЛЬНОЇ ЗАДАЧІ ---")
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
from scipy import stats

data = [
    21, 17, 27, 20, 22, 12, 24, 13, 20, 19, 22, 16, 22,
    9, 21, 16, 23, 16, 21, 24, 18, 11, 22, 15, 23, 21,
    10, 15, 18, 15, 21, 14, 15, 18, 22, 15, 17, 19, 17,
    18, 17, 24, 18, 19, 16, 17, 15, 17, 25, 16
]

x = np.array(data)
n = len(x)
print(f"Обсяг вибірки n = {n}")

unique_elements, counts_elements = np.unique(x, return_counts=True)

var_series = pd.DataFrame({
    'Варіанта (xi)': unique_elements,
    'Частота (ni)': counts_elements
})

print("\n1. Варіаційний ряд:")
print(var_series.to_string(index=False))

mean_val = np.mean(x)

median_val = np.median(x)

mode_result = stats.mode(x, keepdims=True)
mode_val = mode_result.mode[0]
mode_count = mode_result.count[0]

std_dev = np.std(x, ddof=1)
cv = (std_dev / mean_val) * 100

print("\n2. Статистичні показники:")
print(f"{'Середня арифметична (x̄)':<35} {mean_val:.4f}")
print(f"{'Медіана (Me)':<35} {median_val:.4f}")
print(f"{'Мода (Mo)':<35} {mode_val} (зустрічається {mode_count} разів)")
print(f"{'Середнє кв. відхилення (σ або S)':<35} {std_dev:.4f}")
print(f"{'Коефіцієнт варіації (V)':<35} {cv:.2f} %")

homogeneity = "ОДНОРІДНА" if cv <= 30 else "НЕОДНОРІДНА"
print(f"-> Оскільки V = {cv:.2f}%, сукупність є {homogeneity} (поріг 30%).")

plt.figure(figsize=(10, 6))

plt.plot(unique_elements, counts_elements, marker='o', linestyle='-', color='b', linewidth=2, markersize=8)

for i, txt in enumerate(counts_elements):
    plt.annotate(txt, (unique_elements[i], counts_elements[i]), 
                 textcoords="offset points", xytext=(0,10), ha='center')

plt.title('Лінійний графік варіаційного ряду (Полігон розподілу)')
plt.xlabel('Кількість зерен у колосі ($x_i$)')
plt.ylabel('Частота ($n_i$)')
plt.xticks(unique_elements)
plt.grid(True, linestyle='--', alpha=0.7)

print("\n3. Графік побудовано.")
plt.show()