import json
from pathlib import Path

import numpy as np
import pandas as pd

ROOT = Path(__file__).parent
d = pd.read_json(ROOT / 'analysis_data.json')
base = ['封面点击率', '点赞', '评论', '收藏', '涨粉', '分享', '人均观看时长']
rates = []
for metric in ['点赞', '评论', '收藏', '涨粉', '分享']:
    name = metric + '/千次观看'
    d[name] = d[metric] / d['观看量'] * 1000
    rates.append(name)
features = base + rates
rng = np.random.default_rng(20261017)
results = []
B = 200000
for outcome in ['观看量', '曝光']:
    X = d[features].rank().to_numpy()
    X -= X.mean(axis=0)
    X /= np.sqrt(np.sum(X**2, axis=0))
    y = d[outcome].rank().to_numpy()
    y -= y.mean()
    y /= np.sqrt(np.sum(y**2))
    observed = y @ X
    hits = np.zeros(len(features), dtype=int)
    for start in range(0, B, 5000):
        batch = min(5000, B-start)
        shuffled = y[np.argsort(rng.random((batch, len(y))), axis=1)]
        null = shuffled @ X
        hits += (np.abs(null) >= np.abs(observed) - 1e-12).sum(axis=0)
    trim = d[d.source_row != 9]
    for j, metric in enumerate(features):
        results.append({'outcome': outcome, 'feature': metric,
                        'spearman': float(observed[j]),
                        'without_top': float(trim[metric].rank().corr(trim[outcome].rank())),
                        'p_permutation': float((hits[j]+1)/(B+1))})
# One family of 24 exploratory tests, BH-FDR correction.
ps = np.array([row['p_permutation'] for row in results])
order = np.argsort(ps)
q = np.empty_like(ps)
running = 1.0
for i in range(len(ps)-1, -1, -1):
    idx = order[i]
    running = min(running, ps[idx]*len(ps)/(i+1))
    q[idx] = running
for i, row in enumerate(results):
    row['q_bh24'] = float(q[i])

# Bootstrap pair resampling. Re-rank within every resample, retaining ties.
for row in results:
    if row['outcome'] != '观看量' or row['feature'] not in ['封面点击率', '评论/千次观看', '分享/千次观看']:
        continue
    n = len(d)
    samples = rng.integers(0, n, (20000, n))
    a = d[row['feature']].to_numpy()[samples]
    b = d[row['outcome']].to_numpy()[samples]
    ar = pd.DataFrame(a).rank(axis=1).to_numpy()
    br = pd.DataFrame(b).rank(axis=1).to_numpy()
    ar -= ar.mean(axis=1, keepdims=True)
    br -= br.mean(axis=1, keepdims=True)
    denom = np.sqrt((ar**2).sum(axis=1)*(br**2).sum(axis=1))
    valid = denom > 0
    boot = (ar[valid]*br[valid]).sum(axis=1)/denom[valid]
    row['bootstrap_ci95'] = np.quantile(boot,[.025,.975]).tolist()

(ROOT / 'correlation_results.json').write_text(json.dumps(results, ensure_ascii=False, indent=2))
print(json.dumps(results, ensure_ascii=False, indent=2))