API tools faq
paste
mirosh111000

Прикладна_економетрика_ЛР№10_Мірошниченко

mirosh111000
Dec 13th, 2025
85
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Python 20.11 KB | None | 0 0
raw download clone embed print report
  1. import math
  2.  
  3. ALPHA = 0.05
  4. M = 4
  5. TRAIN_START, TRAIN_END = 1, 17
  6. TEST_START, TEST_END = 18, 20
  7. FORECAST_START, FORECAST_END = 21, 24
  8.  
  9. T_CRIT_DF12_A005 = 2.178812829667228
  10. F_CRIT_4_12_A005 = 3.259166726901247
  11.  
  12.  
  13. DATA = [
  14.     (1,  52.0, 72.0, 13.0, 2.7,  5.0),
  15.     (2,  53.0, 74.0, 12.5, 2.8,  5.5),
  16.     (3,  50.0, 72.0, 12.0, 3.0,  5.0),
  17.     (4,  51.0, 73.0, 11.0, 3.2,  6.0),
  18.     (5,  54.0, 70.0, 10.1, 3.2,  7.0),
  19.     (6,  55.0, 67.0,  9.0, 3.3,  8.0),
  20.     (7,  57.0, 67.0,  8.5, 3.4, 10.0),
  21.     (8,  52.0, 62.0,  8.2, 3.6, 10.0),
  22.     (9,  60.0, 72.0,  8.0, 3.7, 10.5),
  23.     (10, 60.0, 72.0,  5.5, 3.7, 11.0),
  24.     (11, 62.0, 74.0,  5.0, 3.4, 13.0),
  25.     (12, 64.0, 75.0,  4.7, 4.0, 10.0),
  26.     (13, 65.0, 76.0,  4.6, 4.2, 12.0),
  27.     (14, 67.0, 80.0,  4.0, 4.3, 13.0),
  28.     (15, 67.0, 82.0,  4.1, 4.7, 14.0),
  29.     (16, 62.0, 84.0,  4.2, 4.8, 14.5),
  30.     (17, 63.0, 84.0,  4.5, 4.8, 15.5),
  31.     (18, 66.0, 87.0,  4.0, 4.9, 17.0),
  32.     (19, 68.0, 90.0,  4.0, 5.0, 16.5),
  33.     (20, 70.0, 92.0,  3.0, 4.7, 17.5),
  34.     (21, None, 92.0,  4.0, 5.2, 17.6),
  35.     (22, None, 93.0,  5.0, 5.3, 17.7),
  36.     (23, None, 93.0,  5.0, 5.4, 17.8),
  37.     (24, None, 94.0,  6.0, 5.4, 17.9),
  38. ]
  39.  
  40.  
  41. def fmt(x, nd=6):
  42.     if x is None:
  43.         return "None"
  44.     if isinstance(x, (int,)) and nd == 0:
  45.         return str(x)
  46.     if isinstance(x, float):
  47.         return f"{x:.{nd}f}"
  48.     return str(x)
  49.  
  50.  
  51. def line(title, w=60):
  52.     print("\n" + title)
  53.     print("-" * w)
  54.  
  55.  
  56. def print_table(headers, rows, nd=6):
  57.     srows = []
  58.     for r in rows:
  59.         srows.append([fmt(v, nd) if isinstance(v, float) else fmt(v) for v in r])
  60.  
  61.     widths = [len(h) for h in headers]
  62.     for r in srows:
  63.         for j, cell in enumerate(r):
  64.             widths[j] = max(widths[j], len(cell))
  65.  
  66.     def row_to_str(rr):
  67.         return " | ".join(rr[j].rjust(widths[j]) for j in range(len(headers)))
  68.  
  69.     print(row_to_str(headers))
  70.     print("-+-".join("-" * widths[j] for j in range(len(headers))))
  71.     for r in srows:
  72.         print(row_to_str(r))
  73.  
  74.  
  75. def mat_transpose(A):
  76.     return [list(row) for row in zip(*A)]
  77.  
  78.  
  79. def mat_mul(A, B):
  80.     r, k = len(A), len(A[0])
  81.     k2, c = len(B), len(B[0])
  82.     if k != k2:
  83.         raise ValueError("mat_mul: несумісні розміри")
  84.     C = [[0.0] * c for _ in range(r)]
  85.     for i in range(r):
  86.         for j in range(c):
  87.             s = 0.0
  88.             for t in range(k):
  89.                 s += A[i][t] * B[t][j]
  90.             C[i][j] = s
  91.     return C
  92.  
  93.  
  94. def mat_vec_mul(A, v):
  95.     r, c = len(A), len(A[0])
  96.     if len(v) != c:
  97.         raise ValueError("mat_vec_mul: несумісні розміри")
  98.     out = [0.0] * r
  99.     for i in range(r):
  100.         s = 0.0
  101.         for j in range(c):
  102.             s += A[i][j] * v[j]
  103.         out[i] = s
  104.     return out
  105.  
  106.  
  107. def identity(n):
  108.     I = [[0.0] * n for _ in range(n)]
  109.     for i in range(n):
  110.         I[i][i] = 1.0
  111.     return I
  112.  
  113.  
  114. def mat_inverse(A):
  115.     n = len(A)
  116.     if n == 0 or len(A[0]) != n:
  117.         raise ValueError("mat_inverse: матриця має бути квадратною")
  118.  
  119.     aug = [A[i][:] + identity(n)[i][:] for i in range(n)]
  120.  
  121.     for col in range(n):
  122.         pivot = col
  123.         best = abs(aug[col][col])
  124.         for r in range(col + 1, n):
  125.             val = abs(aug[r][col])
  126.             if val > best:
  127.                 best = val
  128.                 pivot = r
  129.  
  130.         if best < 1e-15:
  131.             raise ValueError("mat_inverse: матриця вироджена або майже вироджена")
  132.  
  133.         if pivot != col:
  134.             aug[col], aug[pivot] = aug[pivot], aug[col]
  135.  
  136.         piv = aug[col][col]
  137.         invp = 1.0 / piv
  138.         for j in range(2 * n):
  139.             aug[col][j] *= invp
  140.  
  141.         for r in range(n):
  142.             if r == col:
  143.                 continue
  144.             factor = aug[r][col]
  145.             if abs(factor) < 1e-18:
  146.                 continue
  147.             for j in range(2 * n):
  148.                 aug[r][j] -= factor * aug[col][j]
  149.  
  150.     inv = [[aug[i][j] for j in range(n, 2 * n)] for i in range(n)]
  151.     return inv
  152.  
  153.  
  154. def to_rows(data):
  155.     rows = []
  156.     for (month, y, x1, x2, x3, x4) in data:
  157.         rows.append({
  158.             "month": month,
  159.             "Y": y,
  160.             "X1": float(x1),
  161.             "X2": float(x2),
  162.             "X3": float(x3),
  163.             "X4": float(x4),
  164.         })
  165.     return rows
  166.  
  167.  
  168. def slice_rows(rows, a, b):
  169.     return [r for r in rows if a <= r["month"] <= b]
  170.  
  171.  
  172. def mean(vals):
  173.     return sum(vals) / len(vals)
  174.  
  175.  
  176. def build_XY_linear(rows):
  177.     X = []
  178.     y = []
  179.     for r in rows:
  180.         if r["Y"] is None:
  181.             continue
  182.         X.append([1.0, r["X1"], r["X2"], r["X3"], r["X4"]])
  183.         y.append(float(r["Y"]))
  184.     return X, y
  185.  
  186.  
  187. def ols_fit(X, y):
  188.     Xt = mat_transpose(X)
  189.     XtX = mat_mul(Xt, X)
  190.     XtX_inv = mat_inverse(XtX)
  191.  
  192.     XtY = []
  193.     for i in range(len(Xt)):
  194.         s = 0.0
  195.         for k in range(len(y)):
  196.             s += Xt[i][k] * y[k]
  197.         XtY.append([s])
  198.  
  199.     beta_col = mat_mul(XtX_inv, XtY)
  200.     beta = [beta_col[i][0] for i in range(len(beta_col))]
  201.     return beta, XtX_inv
  202.  
  203.  
  204. def predict_linear(rows, beta):
  205.     out = []
  206.     for r in rows:
  207.         x = [1.0, r["X1"], r["X2"], r["X3"], r["X4"]]
  208.         yhat = 0.0
  209.         for j in range(len(beta)):
  210.             yhat += beta[j] * x[j]
  211.         out.append(yhat)
  212.     return out
  213.  
  214.  
  215. def regression_stats(y, yhat, m):
  216.     n = len(y)
  217.     ybar = mean(y)
  218.  
  219.     sst = sum((yi - ybar) ** 2 for yi in y)
  220.     sse = sum((y[i] - yhat[i]) ** 2 for i in range(n))
  221.     ssr = sum((yhat[i] - ybar) ** 2 for i in range(n))
  222.  
  223.     r2 = 1.0 - (sse / sst) if sst != 0.0 else 0.0
  224.     R = math.sqrt(max(0.0, r2))
  225.  
  226.     denom = (n - m - 1)
  227.     r2_adj = 1.0 - (1.0 - r2) * (n - 1) / denom if denom > 0 else float("nan")
  228.  
  229.     Su = math.sqrt(sse / denom) if denom > 0 else float("nan")
  230.  
  231.     return {
  232.         "n": n, "m": m,
  233.         "ybar": ybar,
  234.         "SST": sst, "SSE": sse, "SSR": ssr,
  235.         "R": R, "R2": r2, "R2_adj": r2_adj,
  236.         "Su": Su
  237.     }
  238.  
  239.  
  240. def anova_and_F(stats):
  241.     n = stats["n"]
  242.     m = stats["m"]
  243.     ssr = stats["SSR"]
  244.     sse = stats["SSE"]
  245.  
  246.     df_reg = m
  247.     df_res = n - m - 1
  248.     df_tot = n - 1
  249.  
  250.     msr = ssr / df_reg if df_reg > 0 else float("nan")
  251.     mse = sse / df_res if df_res > 0 else float("nan")
  252.     F = msr / mse if mse != 0.0 else float("inf")
  253.  
  254.     return {
  255.         "df_reg": df_reg, "df_res": df_res, "df_tot": df_tot,
  256.         "msr": msr, "mse": mse, "F": F
  257.     }
  258.  
  259.  
  260. def t_tests(beta, XtX_inv, Su, df_res):
  261.     se = []
  262.     tvals = []
  263.     for i in range(len(beta)):
  264.         v = Su * math.sqrt(max(0.0, XtX_inv[i][i]))
  265.         se.append(v)
  266.         tvals.append(beta[i] / v if v != 0.0 else float("inf"))
  267.     return se, tvals
  268.  
  269.  
  270. def build_XY_power(rows):
  271.     X = []
  272.     y = []
  273.     for r in rows:
  274.         if r["Y"] is None:
  275.             continue
  276.         X.append([1.0, math.log(r["X1"]), math.log(r["X2"]), math.log(r["X3"]), math.log(r["X4"])])
  277.         y.append(math.log(float(r["Y"])))
  278.     return X, y
  279.  
  280.  
  281. def predict_power(rows, beta_log):
  282.     out = []
  283.     for r in rows:
  284.         lnY = beta_log[0]
  285.         lnY += beta_log[1] * math.log(r["X1"])
  286.         lnY += beta_log[2] * math.log(r["X2"])
  287.         lnY += beta_log[3] * math.log(r["X3"])
  288.         lnY += beta_log[4] * math.log(r["X4"])
  289.         out.append(math.exp(lnY))
  290.     return out
  291.  
  292.  
  293. def econ_linear(train_rows, beta):
  294.     y = [r["Y"] for r in train_rows]
  295.     ybar = mean(y)
  296.  
  297.     xbar = {
  298.         "X1": mean([r["X1"] for r in train_rows]),
  299.         "X2": mean([r["X2"] for r in train_rows]),
  300.         "X3": mean([r["X3"] for r in train_rows]),
  301.         "X4": mean([r["X4"] for r in train_rows]),
  302.     }
  303.  
  304.     P = {k: (ybar / xbar[k]) for k in xbar}
  305.  
  306.     marg = {"X1": beta[1], "X2": beta[2], "X3": beta[3], "X4": beta[4]}
  307.  
  308.     E = {k: (marg[k] * xbar[k] / ybar) for k in xbar}
  309.  
  310.     A = sum(E[k] for k in ["X1", "X2", "X3", "X4"])
  311.  
  312.     return ybar, xbar, P, marg, E, A
  313.  
  314.  
  315. def econ_power(train_rows, beta_log):
  316.     y = [r["Y"] for r in train_rows]
  317.     ybar = mean(y)
  318.  
  319.     xbar = {
  320.         "X1": mean([r["X1"] for r in train_rows]),
  321.         "X2": mean([r["X2"] for r in train_rows]),
  322.         "X3": mean([r["X3"] for r in train_rows]),
  323.         "X4": mean([r["X4"] for r in train_rows]),
  324.     }
  325.  
  326.     E = {"X1": beta_log[1], "X2": beta_log[2], "X3": beta_log[3], "X4": beta_log[4]}
  327.  
  328.     marg = {k: (E[k] * ybar / xbar[k]) for k in xbar}
  329.  
  330.     P = {k: (ybar / xbar[k]) for k in xbar}
  331.  
  332.     A = E["X1"] + E["X2"] + E["X3"] + E["X4"]
  333.  
  334.     return ybar, xbar, P, marg, E, A
  335.  
  336.  
  337. def substitution_linear(beta):
  338.     a = {"X1": beta[1], "X2": beta[2], "X3": beta[3], "X4": beta[4]}
  339.     keys = ["X1", "X2", "X3", "X4"]
  340.     hk = {}
  341.     for k in keys:
  342.         for j in keys:
  343.             if k == j:
  344.                 continue
  345.             denom = a[k]
  346.             hk[f"h{k[-1]}{j[-1]}"] = (-a[j] / denom) if denom != 0.0 else float("nan")
  347.     return hk
  348.  
  349.  
  350. def substitution_power(beta_log, train_rows):
  351.     E = {"X1": beta_log[1], "X2": beta_log[2], "X3": beta_log[3], "X4": beta_log[4]}
  352.     xbar = {
  353.         "X1": mean([r["X1"] for r in train_rows]),
  354.         "X2": mean([r["X2"] for r in train_rows]),
  355.         "X3": mean([r["X3"] for r in train_rows]),
  356.         "X4": mean([r["X4"] for r in train_rows]),
  357.     }
  358.     keys = ["X1", "X2", "X3", "X4"]
  359.     hk = {}
  360.     for k in keys:
  361.         for j in keys:
  362.             if k == j:
  363.                 continue
  364.             denom = E[k]
  365.             hk[f"h{k[-1]}{j[-1]}"] = (-(E[j] / denom) * (xbar[k] / xbar[j])) if denom != 0.0 else float("nan")
  366.     return hk
  367.  
  368.  
  369. def forecast_metrics(y_true, y_pred):
  370.     n = len(y_true)
  371.     err = [y_pred[i] - y_true[i] for i in range(n)]
  372.     mae = sum(abs(e) for e in err) / n
  373.     mse = sum(e * e for e in err) / n
  374.     mape = (sum(abs(err[i]) / y_true[i] for i in range(n)) / n) * 100.0
  375.  
  376.     rmse = math.sqrt(sum((y_pred[i] - y_true[i]) ** 2 for i in range(n)) / n)
  377.     denom = math.sqrt(sum(y_true[i] ** 2 for i in range(n)) / n) + math.sqrt(sum(y_pred[i] ** 2 for i in range(n)) / n)
  378.     theil = rmse / denom if denom != 0.0 else float("nan")
  379.  
  380.     return mae, mse, mape, theil
  381.  
  382.  
  383. def main():
  384.     rows_all = to_rows(DATA)
  385.  
  386.     print("=== ЕКОНОМЕТРИЧНА МОДЕЛЬ ПРОДУКТИВНОСТІ ПРАЦІ ===")
  387.     print("Y  – продуктивність праці")
  388.     print("X1 – фондомісткість продукції")
  389.     print("X2 – коефіцієнт плинності робочої сили")
  390.     print("X3 – процент втрат робочого часу")
  391.     print("X4 – стаж роботи")
  392.  
  393.     train_rows = slice_rows(rows_all, TRAIN_START, TRAIN_END)
  394.     test_rows = slice_rows(rows_all, TEST_START, TEST_END)
  395.     future_rows = slice_rows(rows_all, FORECAST_START, FORECAST_END)
  396.  
  397.     line("Вихідні дані (перші/останні рядки)")
  398.     head = rows_all[:6]
  399.     tail = rows_all[-6:]
  400.     out = []
  401.     for r in head + tail:
  402.         out.append([r["month"], r["Y"], r["X1"], r["X2"], r["X3"], r["X4"]])
  403.     print_table(["month", "Y", "X1", "X2", "X3", "X4"], out, nd=3)
  404.  
  405.     line("Лінійна модель (МНК) на 1–17")
  406.     X_lin, y_lin = build_XY_linear(train_rows)
  407.     beta_lin, XtX_inv_lin = ols_fit(X_lin, y_lin)
  408.     yhat_lin_train = mat_vec_mul(X_lin, beta_lin)
  409.  
  410.     stats_lin = regression_stats(y_lin, yhat_lin_train, M)
  411.     anova_lin = anova_and_F(stats_lin)
  412.     se_lin, t_lin = t_tests(beta_lin, XtX_inv_lin, stats_lin["Su"], anova_lin["df_res"])
  413.  
  414.     print("Рівняння:")
  415.     print(f"Ŷ = {fmt(beta_lin[0], 4)} + {fmt(beta_lin[1], 4)}·X1 {fmt(beta_lin[2], 4)}·X2 {fmt(beta_lin[3], 4)}·X3 {fmt(beta_lin[4], 4)}·X4")
  416.  
  417.     print("\nРегресійна статистика:")
  418.     print_table(
  419.         ["Показник", "Значення"],
  420.         [
  421.             ["R", stats_lin["R"]],
  422.             ["R^2", stats_lin["R2"]],
  423.             ["R^2 (Амемія)", stats_lin["R2_adj"]],
  424.             ["n", stats_lin["n"]],
  425.             ["m", stats_lin["m"]],
  426.             ["S_u", stats_lin["Su"]],
  427.         ],
  428.         nd=9
  429.     )
  430.  
  431.     print("\nДисперсійний аналіз (ANOVA):")
  432.     print_table(
  433.         ["Джерело", "df", "SS", "MS"],
  434.         [
  435.             ["Регресії", anova_lin["df_reg"], stats_lin["SSR"], anova_lin["msr"]],
  436.             ["Залишків", anova_lin["df_res"], stats_lin["SSE"], anova_lin["mse"]],
  437.             ["Загальна", anova_lin["df_tot"], stats_lin["SST"], ""],
  438.         ],
  439.         nd=9
  440.     )
  441.  
  442.     print("\nПеревірка значущості моделі (F-тест):")
  443.     print(f"F = {fmt(anova_lin['F'], 6)};  F_кр(α=0.05, df1=4, df2=12) = {fmt(F_CRIT_4_12_A005, 6)}")
  444.     print("Висновок:", "модель значуща" if anova_lin["F"] > F_CRIT_4_12_A005 else "модель незначуща")
  445.  
  446.     print("\nЗначущість параметрів (t-тест, df=12):")
  447.     print(f"t_кр(α=0.05, df=12) = {fmt(T_CRIT_DF12_A005, 6)}")
  448.     trows = []
  449.     for i, name in enumerate(["a0", "a1", "a2", "a3", "a4"]):
  450.         sig = "значущий" if abs(t_lin[i]) > T_CRIT_DF12_A005 else "незначущий"
  451.         trows.append([name, beta_lin[i], se_lin[i], t_lin[i], sig])
  452.     print_table(["парам", "оцінка", "S(a)", "t", "висновок"], trows, nd=9)
  453.  
  454.     line("Степенева модель через логарифмування (МНК) на 1–17")
  455.     X_log, y_log = build_XY_power(train_rows)
  456.     beta_log, XtX_inv_log = ols_fit(X_log, y_log)
  457.  
  458.     yhat_log_train = mat_vec_mul(X_log, beta_log)
  459.     stats_log = regression_stats(y_log, yhat_log_train, M)
  460.     anova_log = anova_and_F(stats_log)
  461.     se_log, t_log = t_tests(beta_log, XtX_inv_log, stats_log["Su"], anova_log["df_res"])
  462.  
  463.     a0_pow = math.exp(beta_log[0])
  464.  
  465.     print("Рівняння (у логарифмах):")
  466.     print(f"lnŶ = {fmt(beta_log[0], 6)} + {fmt(beta_log[1], 6)}·lnX1 + {fmt(beta_log[2], 6)}·lnX2 + {fmt(beta_log[3], 6)}·lnX3 + {fmt(beta_log[4], 6)}·lnX4")
  467.     print("Еквівалент у вихідному масштабі:")
  468.     print(f"Ŷ = {fmt(a0_pow, 6)} · X1^{fmt(beta_log[1], 6)} · X2^{fmt(beta_log[2], 6)} · X3^{fmt(beta_log[3], 6)} · X4^{fmt(beta_log[4], 6)}")
  469.  
  470.     print("\nРегресійна статистика (для lnY):")
  471.     print_table(
  472.         ["Показник", "Значення"],
  473.         [
  474.             ["R", stats_log["R"]],
  475.             ["R^2", stats_log["R2"]],
  476.             ["R^2 (Амемія)", stats_log["R2_adj"]],
  477.             ["n", stats_log["n"]],
  478.             ["m", stats_log["m"]],
  479.             ["S_u", stats_log["Su"]],
  480.         ],
  481.         nd=9
  482.     )
  483.  
  484.     print("\nПеревірка значущості моделі (F-тест):")
  485.     print(f"F = {fmt(anova_log['F'], 6)};  F_кр(α=0.05, df1=4, df2=12) = {fmt(F_CRIT_4_12_A005, 6)}")
  486.     print("Висновок:", "модель значуща" if anova_log["F"] > F_CRIT_4_12_A005 else "модель незначуща")
  487.  
  488.     print("\nЗначущість параметрів (t-тест, df=12):")
  489.     print(f"t_кр(α=0.05, df=12) = {fmt(T_CRIT_DF12_A005, 6)}")
  490.     trows = []
  491.     for i, name in enumerate(["b0", "b1", "b2", "b3", "b4"]):
  492.         sig = "значущий" if abs(t_log[i]) > T_CRIT_DF12_A005 else "незначущий"
  493.         trows.append([name, beta_log[i], se_log[i], t_log[i], sig])
  494.     print_table(["парам", "оцінка", "S(b)", "t", "висновок"], trows, nd=9)
  495.  
  496.     line("Економічні характеристики (на середніх значеннях 1–17)")
  497.  
  498.     ybar_l, xbar_l, P_l, marg_l, E_l, A_l = econ_linear(train_rows, beta_lin)
  499.     print("Лінійна модель:")
  500.     rows_e = []
  501.     for k in ["X1", "X2", "X3", "X4"]:
  502.         rows_e.append([k, xbar_l[k], P_l[k], marg_l[k], E_l[k]])
  503.     print_table(["чинник", "X̄", "Pj=Ȳ/X̄", "dY/dX", "E (на середніх)"], rows_e, nd=9)
  504.     print(f"Загальна еластичність A = ΣE = {fmt(A_l, 6)}")
  505.     print(f"Інтерпретація: якщо всі чинники зростають на 1%, Y зміниться приблизно на {fmt(A_l, 4)}% (для лінійної).")
  506.  
  507.     sub_l = substitution_linear(beta_lin)
  508.     print("\nНорми заміщення (лінійна): hk_j = -a_j/a_k")
  509.     keys_order = ["h12", "h13", "h14", "h21", "h23", "h24", "h31", "h32", "h34", "h41", "h42", "h43"]
  510.     sub_rows = [[k, sub_l[k]] for k in keys_order]
  511.     print_table(["показник", "значення"], sub_rows, nd=6)
  512.  
  513.     ybar_p, xbar_p, P_p, marg_p, E_p, A_p = econ_power(train_rows, beta_log)
  514.     print("\nСтепенева модель:")
  515.     rows_e = []
  516.     for k in ["X1", "X2", "X3", "X4"]:
  517.         rows_e.append([k, xbar_p[k], P_p[k], marg_p[k], E_p[k]])
  518.     print_table(["чинник", "X̄", "Pj=Ȳ/X̄", "dY/dX (на середніх)", "E=a_j"], rows_e, nd=9)
  519.     print(f"Загальна еластичність A = Σa_j = {fmt(A_p, 6)}")
  520.     print(f"Інтерпретація: якщо всі чинники зростають на 1%, Y зміниться приблизно на {fmt(A_p, 4)}% (для степеневої).")
  521.  
  522.     sub_p = substitution_power(beta_log, train_rows)
  523.     print("\nНорми заміщення (степенева, на середніх): hk_j = -(a_j/a_k) * (X̄k/X̄j)")
  524.     sub_rows = [[k, sub_p[k]] for k in keys_order]
  525.     print_table(["показник", "значення"], sub_rows, nd=6)
  526.  
  527.     line("Перевірка прогнозу на 18–20 (MAE, MSE, MAPE, Тейл)")
  528.  
  529.     y_test = [r["Y"] for r in test_rows]
  530.     yhat_test_lin = predict_linear(test_rows, beta_lin)
  531.     yhat_test_pow = predict_power(test_rows, beta_log)
  532.  
  533.     cmp_rows = []
  534.     for i in range(len(test_rows)):
  535.         r = test_rows[i]
  536.         cmp_rows.append([
  537.             r["month"],
  538.             y_test[i],
  539.             yhat_test_lin[i],
  540.             (yhat_test_lin[i] - y_test[i]),
  541.             yhat_test_pow[i],
  542.             (yhat_test_pow[i] - y_test[i]),
  543.         ])
  544.     print_table(["month", "Y", "Ŷ_lin", "u_lin", "Ŷ_pow", "u_pow"], cmp_rows, nd=6)
  545.  
  546.     mae_l, mse_l, mape_l, theil_l = forecast_metrics(y_test, yhat_test_lin)
  547.     mae_p, mse_p, mape_p, theil_p = forecast_metrics(y_test, yhat_test_pow)
  548.  
  549.     print("\nПоказники якості прогнозу:")
  550.     print_table(
  551.         ["Модель", "MAE", "MSE", "MAPE,%", "Тейл"],
  552.         [
  553.             ["Лінійна", mae_l, mse_l, mape_l, theil_l],
  554.             ["Степенева", mae_p, mse_p, mape_p, theil_p],
  555.         ],
  556.         nd=9
  557.     )
  558.  
  559.     line("Прогноз Ŷ на 21–24 (за очікуваними X)")
  560.  
  561.     yhat_future_lin = predict_linear(future_rows, beta_lin)
  562.     yhat_future_pow = predict_power(future_rows, beta_log)
  563.  
  564.     out_rows = []
  565.     for i in range(len(future_rows)):
  566.         r = future_rows[i]
  567.         out_rows.append([
  568.             r["month"], r["X1"], r["X2"], r["X3"], r["X4"],
  569.             yhat_future_lin[i], yhat_future_pow[i], (yhat_future_pow[i] - yhat_future_lin[i])
  570.         ])
  571.     print_table(["month", "X1", "X2", "X3", "X4", "Ŷ_lin", "Ŷ_pow", "Δ(pow-lin)"], out_rows, nd=6)
  572.  
  573.     line("Висновки")
  574.     print("1) Значущість моделей (F):")
  575.     print(f"   Лінійна:   F={fmt(anova_lin['F'], 4)} > Fкр={fmt(F_CRIT_4_12_A005, 4)} -> значуща")
  576.     print(f"   Степенева: F={fmt(anova_log['F'], 4)} > Fкр={fmt(F_CRIT_4_12_A005, 4)} -> значуща")
  577.  
  578.     print("\n2) Значущі параметри (t, α=0.05, df=12):")
  579.     sig_lin = [("a"+str(i), abs(t_lin[i]) > T_CRIT_DF12_A005) for i in range(5)]
  580.     sig_pow = [("b"+str(i), abs(t_log[i]) > T_CRIT_DF12_A005) for i in range(5)]
  581.     print("   Лінійна:", ", ".join([f"{n}:{'+' if ok else '-'}" for n, ok in sig_lin]))
  582.     print("   Степенева:", ", ".join([f"{n}:{'+' if ok else '-'}" for n, ok in sig_pow]))
  583.  
  584.     print("\n3) Загальна еластичність:")
  585.     print(f"   Лінійна A = {fmt(A_l, 6)}  (≈ {fmt(A_l, 4)}% при +1% усіх X)")
  586.     print(f"   Степенева A = {fmt(A_p, 6)} (≈ {fmt(A_p, 4)}% при +1% усіх X)")
  587.  
  588.  
  589.  
  590. if __name__ == "__main__":
  591.     main()
  592.  
Advertisement
Add Comment
Please, Sign In to add comment
Public Pastes
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!