# Data sources:
# New study permits: 40k-60k
# New work permits: 280k-360k
# Study permits (new): capped at 437k for 2025. H1 likely 40k-60k new approvals (cap + seasonality).
# Work permits (new): ~825k finalized Jan-Jul, but 40-50% are renewals. H1 ~280k-360k new approvals.
# Combined H1 2025 “new” study+work ≈ 320k-420k.
# StatsCan shows NPR stock fell in Q1 2025, meaning exits/extensions outweighed new arrivals.
# https://www.canada.ca/en/immigration-refugees-citizenship/corporate/mandate/corporate-initiatives/levels/inventories-backlogs.html
# https://www.canada.ca/en/immigration-refugees-citizenship/news/notices/2025-provincial-territorial-allocations-under-international-student-cap.html https://open.canada.ca/data/en/dataset/9b34e712-513f-44e9-babf-9df4f7256550
# https://open.canada.ca/data/en/dataset/90115b00-f9b8-49e8-afa3-b4cff8facaee
# https://open.canada.ca/data/en/dataset/360024f2-17e9-4558-bfc1-3616485d65b9
# https://www150.statcan.gc.ca/t1/tbl1/en/tv.action?pid=1710012101
# https://www150.statcan.gc.ca/n1/daily-quotidien/250618/dq250618a-eng.htm
# Births: https://www150.statcan.gc.ca/t1/tbl1/en/dtl!downloadDbLoadingData-nonTraduit.action?pid=1310041501&latestN=0&startDate=20000101&endDate=99990101&csvLocale=en&selectedMembers=%5B%5B1%5D%2C%5B%5D%2C%5B1%5D%5D&checkedLevels=1D1%2C1D2 (https://www150.statcan.gc.ca/t1/tbl1/en/tv.action?pid=1310041501) - 351,878
# Deaths: https://www150.statcan.gc.ca/t1/tbl1/en/dtl!downloadDbLoadingData-nonTraduit.action?pid=1310070801&latestN=0&startDate=20000101&endDate=99990101&csvLocale=en&selectedMembers=%5B%5B1%5D%2C%5B1%2C2%2C3%2C4%2C5%2C6%2C7%2C8%2C9%2C10%2C11%2C12%2C13%5D%2C%5B1%2C2%5D%5D&checkedLevels= (https://www150.statcan.gc.ca/t1/tbl1/en/tv.action?pid=1310070801) - 326,215

import matplotlib.pyplot as plt
from matplotlib.lines import Line2D

TOTAL_POP = 41_290_000
BIRTHS = 351_878
DEATHS = 326_215
OBSERVED_POP_CHANGE = 20_107
NATURAL_INCREASE = BIRTHS - DEATHS

NEW_STUDY_RANGE = (40_000, 60_000)
NEW_WORK_RANGE  = (280_000, 360_000)
RENEWED_WORK_RANGE = (330_000, 410_000)

STUDENT_EXPIRY_SHARE_RANGE = (0.30, 0.40)
WORK_EXPIRY_SHARE_RANGE    = (0.60, 0.70)

def midpoint(lo: int, hi: int) -> float:
    return (lo + hi) / 2.0

def halfspan(lo: int, hi: int) -> float:
    return (hi - lo) / 2.0

NEW_VISAS_RANGE = (NEW_STUDY_RANGE[0] + NEW_WORK_RANGE[0],
                   NEW_STUDY_RANGE[1] + NEW_WORK_RANGE[1])

EXPIRIES_TOTAL_RANGE = (
    NATURAL_INCREASE + NEW_VISAS_RANGE[0] + RENEWED_WORK_RANGE[0] - OBSERVED_POP_CHANGE,
    NATURAL_INCREASE + NEW_VISAS_RANGE[1] + RENEWED_WORK_RANGE[1] - OBSERVED_POP_CHANGE
)

EXPIRING_STUDENTS_RANGE = (
    int(EXPIRIES_TOTAL_RANGE[0] * STUDENT_EXPIRY_SHARE_RANGE[0]),
    int(EXPIRIES_TOTAL_RANGE[1] * STUDENT_EXPIRY_SHARE_RANGE[1]),
)
EXPIRING_WORK_RANGE = (
    int(EXPIRIES_TOTAL_RANGE[0] * WORK_EXPIRY_SHARE_RANGE[0]),
    int(EXPIRIES_TOTAL_RANGE[1] * WORK_EXPIRY_SHARE_RANGE[1]),
)

categories = [
    "Births",
    "Deaths",
    "New Study Visas",
    "New Work Visas",
    "Renewed Work Permits",
    "Expiring Student Visas",
    "Expiring Work Permits",
    "Population Change",
]

values = [
    float(BIRTHS),
    -float(DEATHS),
    midpoint(*NEW_STUDY_RANGE),
    midpoint(*NEW_WORK_RANGE),
    midpoint(*RENEWED_WORK_RANGE),
    -midpoint(*EXPIRING_STUDENTS_RANGE),
    -midpoint(*EXPIRING_WORK_RANGE),
    float(OBSERVED_POP_CHANGE),
]

yerr = [
    0.0,
    0.0,
    halfspan(*NEW_STUDY_RANGE),
    halfspan(*NEW_WORK_RANGE),
    halfspan(*RENEWED_WORK_RANGE),
    halfspan(*EXPIRING_STUDENTS_RANGE),
    halfspan(*EXPIRING_WORK_RANGE),
    0.0,
]

colors = [
    "#2ca02c",  # Births
    "#d62728",  # Deaths
    "#17becf",  # New Study (teal)
    "#1f77b4",  # New Work (blue)
    "#9467bd",  # Renewed Work (purple)
    "#ff7f0e",  # Expiring Student (orange)
    "#8c564b",  # Expiring Work (brown)
    "#7f7f7f",  # Pop change (grey)
]

plt.figure(figsize=(14, 7.5))
bars = plt.bar(categories, values, color=colors)

for i, (v, e) in enumerate(zip(values, yerr)):
    if e > 0:
        plt.errorbar(i, v, yerr=e, fmt='none', capsize=6, elinewidth=1.5, ecolor='black')
        plt.scatter([i], [v], s=40, marker='D', edgecolor='black', zorder=3)

pad = 120_000
for bar, v, e in zip(bars, values, yerr):
    x = bar.get_x() + bar.get_width() / 2
    va = 'bottom' if v >= 0 else 'top'
    if e > 0:
        lo, hi = v - e, v + e
        txt = f"central estimate*: {int(round(v)):,}\nrange: {int(round(lo)):,}-{int(round(hi)):,}"
    else:
        txt = f"{int(round(abs(v))):,}"
    y = v + (pad if v >= 0 else -pad)
    plt.text(x, y, txt, ha='center', va=va, fontsize=9)

def add_group_bracket(ax, x0, x1, y, label, invert=False):
    if not invert:
        ax.plot([x0, x0, x1, x1], [y, y+40_000, y+40_000, y], lw=1.5, color='black')
        ax.text((x0+x1)/2, y+55_000, label, ha='center', va='bottom', fontsize=10)
    else:
        ax.plot([x0, x0, x1, x1], [y, y-40_000, y-40_000, y], lw=1.5, color='black')
        ax.text((x0+x1)/2, y-55_000, label, ha='center', va='top', fontsize=10)

ax = plt.gca()

ns_idx, nw_idx, rw_idx = 2, 3, 4
expS_idx, expW_idx = 5, 6

y_ns = values[ns_idx] + yerr[ns_idx]
y_nw = values[nw_idx] + yerr[nw_idx]
y_rw = values[rw_idx] + yerr[rw_idx]
y_expS = values[expS_idx] - yerr[expS_idx]  # negative
y_expW = values[expW_idx] - yerr[expW_idx]

y_group_new = max(y_ns, y_nw) + 180_000
sum_new = int(round(values[ns_idx] + values[nw_idx]))
add_group_bracket(ax, ns_idx, nw_idx, y_group_new, f"New visas (Study+Work sum): {sum_new:,}")

y_group_new_renewed = max(y_group_new + 220_000, max(y_ns, y_nw, y_rw) + 260_000)
sum_new_renewed = int(round(values[ns_idx] + values[nw_idx] + values[rw_idx]))
add_group_bracket(ax, ns_idx, rw_idx, y_group_new_renewed, f"New (Study+Work) + Renewed Work (sum): {sum_new_renewed:,}")

y_group_exp = min(y_expS, y_expW) - 220_000
sum_exp = int(round(abs(values[expS_idx] + values[expW_idx])))
add_group_bracket(ax, expS_idx, expW_idx, y_group_exp, f"All expiries (sum): {sum_exp:,}", invert=True)

legend_elems = [
    Line2D([0], [0], color='black', lw=1.5, label='Range (vertical line)'),
    Line2D([0], [0], marker='D', color='w', markeredgecolor='black', label='Central estimate*', markersize=7),
]
plt.legend(handles=legend_elems, loc='upper right')

plt.axhline(0, linewidth=0.9, color='black')
plt.title("Population Components: Canada (H1 2025)\nNew visas split into Study vs Work; ranges shown with central estimate*", pad=16)
plt.ylabel("Count (people)")

upper = max((v + e) for v, e in zip(values, yerr)) * 1.75 + 200_000
lower = min((v - e) for v, e in zip(values, yerr)) * 1.75 - 250_000
plt.ylim(lower, upper)

plt.xticks(rotation=10)
plt.tight_layout()

plt.figtext(
    0.5, 0.01,
    "*Estimate based on IRCC/StatCan data and a heuristic expiry split (30-40% students, 60-70% workers).",
    ha="center", fontsize=9, style="italic"
)

plt.show()