# -*- coding: utf-8 -*-
"""AI 확산 지수 — 특허·논문·뉴스로 만드는 총량 지수(분기)와 산업별 노출 지수(설문 21개 업종 → 62개 산업).

총량 지수(시계열)
  세 자료원의 로그 수준을 각자 2015–2024년 평균으로 뺀 값(로그 편차)을 분기 단위로 맞추고, 쓸 수 있는 자료원의 평균을 낸다.
    · 특허 : FindIP IPC G06N 5개국 출원(연) → 분기 보간. 2025년은 공개 지연으로 덜 잡혀 쓰지 않는다.
    · 논문 : OpenAlex 'Artificial intelligence' 개념 논문(연) → 분기 보간.
    · 뉴스 : Event Registry 영문 'artificial intelligence' 기사(일) → 분기 합계(2015년부터).
  연간 자료의 분기 보간은 연 중앙(2분기 말)에 로그 수준을 두고 선형보간한 뒤, 해마다 분기 합이 연간 값과 같도록 비례 조정한다.
  AI 확산 갭 = 총량 지수 − HP 추세(λ=1600). 2015년 이전은 0(뉴스 자료 시작 전, AI 확산 이전 시기)으로 둔다.

산업별 노출 지수(횡단)
  · 뉴스 : 2025년 업종 키워드 기사 중 'artificial intelligence'를 함께 다룬 비율(%)
  · 논문 : 2023–2025년 OpenAlex 분야별 AI 논문 비중을 분야→업종 대응표로 가중평균(분야 전체 논문 수 가중)
  · 특허 : 5개국 G06N 특허(2020–2025)의 동시분류를 IPC→업종 대응표로 배분한 건수 ÷ 업종 부가가치(2025년, 조 원).
           동시분류 건수는 상위 20개 IPC 표에, 그 밖의 대응표 IPC를 포함–배제(|G06N| + |X| − |G06N∪X|)로 구해 더한다.
  자료원마다 로그값을 업종 간 표준화(z)하고, 쓸 수 있는 자료원의 평균을 노출 지수로 둔다.
"""
from __future__ import annotations

import numpy as np
import pandas as pd

YEARS = list(range(2010, 2026))


def quarters(y0: int, y1: int) -> list[str]:
    return [f"{y}Q{q}" for y in range(y0, y1 + 1) for q in (1, 2, 3, 4)]


def annual_to_quarterly(s: pd.Series) -> pd.Series:
    """연간 합계 → 분기: 연 중앙에 로그 수준을 두고 분기 간 선형보간, 해마다 분기 합 = 연간 값(비례 조정)."""
    s = s.dropna().astype(float)
    ys = list(s.index)
    qs = quarters(min(ys), max(ys))
    t_ann = np.array([y + 0.5 for y in ys])                     # 연 중앙
    t_q = np.array([int(q[:4]) + (int(q[-1]) - 0.5) / 4 for q in qs])
    lq = np.interp(t_q, t_ann, np.log(s.values / 4.0))
    out = pd.Series(np.exp(lq), index=qs)
    for y in ys:
        idx = [q for q in qs if q.startswith(str(y))]
        out[idx] *= s[y] / out[idx].sum()
    return out


def daily_to_quarterly(d: pd.Series) -> pd.Series:
    q = d.groupby([f"{t.year}Q{(t.month - 1) // 3 + 1}" for t in d.index]).sum()
    return q.sort_index()


def log_dev(s: pd.Series, base_years=(2015, 2024)) -> pd.Series:
    """로그 수준 − 기준기간(분기) 평균 로그 수준."""
    ls = np.log(s.astype(float))
    base = [q for q in ls.index if base_years[0] <= int(q[:4]) <= base_years[1]]
    return ls - ls.loc[base].mean()


def hp_filter(y: np.ndarray, lam: float = 1600.0) -> np.ndarray:
    """호드릭–프레스콧 추세(희소 행렬 없이 밀집 행렬로 계산, 표본이 짧아 충분하다)."""
    n = len(y)
    D = np.zeros((n - 2, n))
    for i in range(n - 2):
        D[i, i:i + 3] = [1.0, -2.0, 1.0]
    return np.linalg.solve(np.eye(n) + lam * D.T @ D, y)


def composite_index(patent_annual: pd.Series, paper_annual: pd.Series, news_daily: pd.Series,
                    patent_last_year: int = 2024) -> pd.DataFrame:
    """세 자료원의 분기 로그 편차와 평균(총량 지수), HP 추세, 갭."""
    pat = log_dev(annual_to_quarterly(patent_annual.loc[:patent_last_year]))
    pap = log_dev(annual_to_quarterly(paper_annual))
    new = log_dev(daily_to_quarterly(news_daily))
    qs = quarters(2015, 2025)
    df = pd.DataFrame({"특허": pat.reindex(qs), "논문": pap.reindex(qs), "뉴스": new.reindex(qs)})
    df["총량지수"] = df[["특허", "논문", "뉴스"]].mean(axis=1, skipna=True)
    df["추세"] = hp_filter(df["총량지수"].values)
    df["갭"] = df["총량지수"] - df["추세"]
    df.index.name = "분기"
    return df


def gap_series(comp: pd.DataFrame, qs: list[str]) -> pd.Series:
    """추정 표본 분기(2000Q1–2025Q4)에 맞춘 AI 확산 갭. 2015년 이전은 0."""
    return comp["갭"].reindex(qs).fillna(0.0).rename("AI확산갭")


# ── 산업별 노출 ─────────────────────────────────────────────────────────
def exposure_news(er: pd.DataFrame) -> pd.Series:
    return er.set_index("설문업종")["AI비중_%"].rename("뉴스")


def exposure_papers(field_share: pd.DataFrame, fmap: pd.DataFrame) -> pd.Series:
    fs = field_share.reset_index().rename(columns={field_share.index.name or "index": "분야(OpenAlex field)"})
    m = fmap.merge(fs, on="분야(OpenAlex field)", how="left")
    m["w"] = m["배분비중"] * m["전체논문"]
    g = m.groupby("설문업종").apply(lambda d: float((d["AI비중_%"] * d["w"]).sum() / d["w"].sum()))
    return g.rename("논문")


def merge_coipc(top20: pd.DataFrame, supplement: pd.DataFrame) -> pd.DataFrame:
    """상위 20개 동시분류 표에 포함–배제로 구한 보완 건수(|G06N|+|X|−|G06N∪X|)를 더한다. 겹치는 IPC는 상위 20개 값을 쓴다."""
    col = top20.columns[1]
    sup = supplement[["IPC", "G06N_동시분류_건수_5개국_2020_2025"]].rename(columns={"G06N_동시분류_건수_5개국_2020_2025": col})
    sup = sup[~sup["IPC"].isin(top20["IPC"])]
    return pd.concat([top20, sup], ignore_index=True)


def exposure_patents(coipc: pd.DataFrame, imap: pd.DataFrame, va_by_ind: pd.Series) -> pd.Series:
    c = coipc.rename(columns={coipc.columns[1]: "건수"})
    m = imap.merge(c, on="IPC", how="left").fillna({"건수": 0})
    m["배분건수"] = m["건수"] * m["배분비중"]
    per = m.groupby("설문업종")["배분건수"].sum()
    return (per / va_by_ind.reindex(per.index)).rename("특허")       # 건 / 조 원


def exposure_index(parts: list[pd.Series], industries: list[str]) -> pd.DataFrame:
    df = pd.concat(parts, axis=1).reindex(industries)
    z = np.log(df.where(df > 0))
    z = (z - z.mean()) / z.std(ddof=0)
    df["노출지수_z"] = z.mean(axis=1, skipna=True)
    df["자료원수"] = z.notna().sum(axis=1)
    return df


MANUF = ["01. 전기전자", "02. 정밀기기", "03. 섬유가죽", "04. 석유화학", "05. 비금속광물",
         "06. 금속제품", "07. 기계장비", "08. 자동차", "09. 식료품"]


def exposure_by_group(z: pd.Series, va_by_group: pd.Series) -> pd.Series:
    """21개 업종 노출 지수에 '전체'·'제조업'을 더한다(대응표에서 설문 업종이 없는 62개 산업용).
    '전체' = 21개 업종의 부가가치 가중평균, '제조업' = 제조 9개 업종의 부가가치 가중평균."""
    w = va_by_group.reindex(z.index).fillna(0.0)
    out = z.copy()
    out["전체"] = float((z * w).sum() / w.sum())
    m = [g for g in MANUF if g in z.index]
    out["제조업"] = float((z[m] * w[m]).sum() / w[m].sum())
    return out


def survey_scale_shape(z62: np.ndarray, a62: np.ndarray, w62: np.ndarray) -> np.ndarray:
    """노출 지수의 순위·간격을 설문 강도의 척도로 옮긴다: a' = ā + s_a·(z − z̄)/s_z (부가가치 가중 평균·표준편차).
    객관 지표는 산업 간 차이가 설문보다 훨씬 커서(정보통신 특허 밀도는 다른 업종의 수십 배) 그대로 쓰면 배분이 한두 산업에 몰린다.
    설문이 말하는 효과의 평균 크기와 산업 간 분산은 유지하고, 어느 산업이 더 노출되었는지만 객관 지표로 바꾸는 방식이다."""
    w = w62 / w62.sum()
    ma, mz = float((w * a62).sum()), float((w * z62).sum())
    sa = float(np.sqrt((w * (a62 - ma) ** 2).sum())); sz = float(np.sqrt((w * (z62 - mz) ** 2).sum()))
    return ma + sa * (z62 - mz) / sz


def spearman(a: pd.Series, b: pd.Series) -> float:
    d = pd.concat([a, b], axis=1).dropna()
    return float(d.rank().corr().iloc[0, 1])