1# Author: Simon-Pierre Boucher — simon-pierre.boucher@uqo.ca2"""Information aggregation and deliberation: closed forms and Monte Carlo (numpy only).34Notation follows `paper/sections/theory.tex`:5 * K agents observe signals s_k = v + eta_k, eta ~ N(0, Sigma).6 * Linear pool w (sum to one): MSE(w) = w' Sigma w.7 * Exchangeable Sigma = sigma^2 [(1-rho) I + rho 11']: MSE_mean = sigma^2 [(1-rho)/K + rho].8 * Opinion-only deliberation: x^{r+1} = (1-lambda) x^r + lambda T(x^r) 1, T = mean or median.9 * Evidence model: J evidence items z_j = v + c + u_j (c common, u_j iid), agent k uses the uniform10 average of its view A_k plus a provider-common error g_{m(k)} and an idiosyncratic reasoning error e_k.11 Evidence exchange enlarges the views; the implied evidence weights omega = (1/K) sum_k a_k determine12 the pooled MSE through sigma_u^2 ||omega||^2.13"""14from __future__ import annotations1516from dataclasses import dataclass, field1718import numpy as np1920# --------------------------------------------------------------------------------------------------21# Static linear pools (Propositions 1–3)22# --------------------------------------------------------------------------------------------------232425def exchangeable_cov(K: int, sigma: float, rho: float) -> np.ndarray:26 """Sigma = sigma^2 [(1-rho) I + rho 11']; requires rho >= -1/(K-1) for positive semi-definiteness."""27 if K > 1 and rho < -1.0 / (K - 1):28 raise ValueError("rho below -1/(K-1): not a valid correlation matrix")29 return sigma**2 * ((1.0 - rho) * np.eye(K) + rho * np.ones((K, K)))303132def mse_mean_closed_form(K: int, sigma: float, rho: float) -> float:33 """Proposition 1: MSE of the equal-weight mean of K exchangeable signals."""34 return sigma**2 * ((1.0 - rho) / K + rho)353637def mse_pool(w: np.ndarray, Sigma: np.ndarray, b: np.ndarray | None = None) -> float:38 """MSE of the linear pool w (weights summing to one) with error covariance Sigma and bias vector b."""39 w = np.asarray(w, dtype=float)40 val = float(w @ Sigma @ w)41 if b is not None:42 val += float(w @ np.asarray(b, dtype=float)) ** 243 return val444546def optimal_weights(Sigma: np.ndarray) -> np.ndarray:47 """Bates–Granger weights w* = Sigma^{-1} 1 / (1' Sigma^{-1} 1) (unbiased signals, known Sigma)."""48 K = Sigma.shape[0]49 x = np.linalg.solve(Sigma, np.ones(K))50 return x / x.sum()515253def mse_optimal(Sigma: np.ndarray) -> float:54 """Proposition 2: 1 / (1' Sigma^{-1} 1)."""55 K = Sigma.shape[0]56 return 1.0 / float(np.ones(K) @ np.linalg.solve(Sigma, np.ones(K)))575859def mse_equal(Sigma: np.ndarray) -> float:60 """1' Sigma 1 / K^2."""61 K = Sigma.shape[0]62 return float(Sigma.sum()) / K**2636465def nested_cov(n_per_provider: int, M: int, sigma_c: float, sigma_g: float, sigma_e: float) -> np.ndarray:66 """Proposition 3: eta_k = c + g_{m(k)} + e_k with K = n*M agents in M provider blocks."""67 K = n_per_provider * M68 block = np.kron(np.eye(M), np.ones((n_per_provider, n_per_provider)))69 return sigma_c**2 * np.ones((K, K)) + sigma_g**2 * block + sigma_e**2 * np.eye(K)707172def mse_nested_closed_form(n_per_provider: int, M: int, sigma_c: float, sigma_g: float, sigma_e: float) -> float:73 K = n_per_provider * M74 return sigma_c**2 + sigma_g**2 / M + sigma_e**2 / K757677def simulate_mean_mse(K: int, sigma: float, rho: float, n_draws: int, rng: np.random.Generator) -> float:78 """Monte Carlo MSE of the equal-weight mean under the exchangeable model (v = 0 without loss)."""79 c = rng.standard_normal(n_draws) * (sigma * np.sqrt(max(rho, 0.0)))80 e = rng.standard_normal((n_draws, K)) * (sigma * np.sqrt(1.0 - rho))81 if rho < 0: # negative rho: draw from the full covariance instead of the factor representation82 L = np.linalg.cholesky(exchangeable_cov(K, sigma, rho))83 eta = rng.standard_normal((n_draws, K)) @ L.T84 else:85 eta = c[:, None] + e86 return float(np.mean(eta.mean(axis=1) ** 2))878889# --------------------------------------------------------------------------------------------------90# Opinion-only deliberation (Propositions 4–5, 7)91# --------------------------------------------------------------------------------------------------929394def deliberate_opinion_only(x0: np.ndarray, rounds: int, lam: float, anchor: str = "median") -> np.ndarray:95 """Apply x^{r+1} = (1-lam) x^r + lam T(x^r) 1 for `rounds` rounds to each row of x0 (draws x K).96 Returns an array (rounds+1, draws, K) with the trajectory."""97 x = np.asarray(x0, dtype=float)98 traj = np.empty((rounds + 1,) + x.shape)99 traj[0] = x100 for r in range(rounds):101 stat = np.median(x, axis=1) if anchor == "median" else x.mean(axis=1)102 x = (1.0 - lam) * x + lam * stat[:, None]103 traj[r + 1] = x104 return traj105106107def degroot_consensus_weights(W: np.ndarray) -> np.ndarray:108 """Left Perron vector pi of a primitive row-stochastic W (pi' W = pi', pi' 1 = 1)."""109 vals, vecs = np.linalg.eig(W.T)110 k = int(np.argmin(np.abs(vals - 1.0)))111 pi = np.real(vecs[:, k])112 return pi / pi.sum()113114115def inverse_variance_weights(Sigma: np.ndarray) -> np.ndarray:116 """Weights proportional to 1/Sigma_kk (ignore correlations), as used by the empirical aggregator."""117 w = 1.0 / np.diag(Sigma)118 return w / w.sum()119120121def deliberation_cov(Sigma: np.ndarray, rounds: int, lam: float) -> np.ndarray:122 """Error covariance after `rounds` rounds of mean-anchored linear deliberation: W^r Sigma W^r."""123 K = Sigma.shape[0]124 W = (1.0 - lam) * np.eye(K) + lam * np.ones((K, K)) / K125 Wr = np.linalg.matrix_power(W, rounds)126 return Wr @ Sigma @ Wr.T127128129# --------------------------------------------------------------------------------------------------130# Evidence model with exchange (Proposition 6) and placebo roles (Proposition 8)131# --------------------------------------------------------------------------------------------------132133134@dataclass135class EvidenceCommittee:136 """K agents, J evidence items, views as index lists; provider blocks; reasoning-error variances."""137138 J: int139 views: list[np.ndarray]140 providers: np.ndarray # provider index m(k)141 tau2: np.ndarray # reasoning-error variance of agent k (prompt-specific), may depend on the slot142 sigma_c: float = 0.0143 sigma_u: float = 1.0144 sigma_g: float = 0.0145146 @property147 def K(self) -> int:148 return len(self.views)149150 def view_matrix(self, views: list[np.ndarray] | None = None) -> np.ndarray:151 """A (K x J) with uniform weights on each view (the BLUE within a view under iid evidence errors)."""152 views = self.views if views is None else views153 A = np.zeros((self.K, self.J))154 for k, idx in enumerate(views):155 A[k, np.asarray(idx)] = 1.0 / len(idx)156 return A157158 def implied_evidence_weights(self, views: list[np.ndarray] | None = None) -> np.ndarray:159 return self.view_matrix(views).mean(axis=0)160161 def mse_equal_closed_form(self, views: list[np.ndarray] | None = None, tau2: np.ndarray | None = None) -> float:162 """sigma_c^2 + sigma_u^2 ||omega||^2 + sigma_g^2 sum_m (n_m/K)^2 + (1/K^2) sum_k tau2_k."""163 omega = self.implied_evidence_weights(views)164 tau2 = self.tau2 if tau2 is None else tau2165 _, counts = np.unique(self.providers, return_counts=True)166 prov = self.sigma_g**2 * float(np.sum((counts / self.K) ** 2))167 return self.sigma_c**2 + self.sigma_u**2 * float(omega @ omega) + prov + float(tau2.sum()) / self.K**2168169 def draw(self, n_draws: int, rng: np.random.Generator, views: list[np.ndarray] | None = None,170 tau2: np.ndarray | None = None, shocks: dict | None = None) -> np.ndarray:171 """Round-0 signals (n_draws x K) with v = 0. `shocks` allows common random numbers across variants."""172 tau2 = self.tau2 if tau2 is None else tau2173 if shocks is None:174 shocks = self.draw_shocks(n_draws, rng)175 A = self.view_matrix(views)176 c, u, g, e = shocks["c"], shocks["u"], shocks["g"], shocks["e"]177 return c[:, None] + u @ A.T + g[:, self.providers] + e * np.sqrt(tau2)[None, :]178179 def draw_shocks(self, n_draws: int, rng: np.random.Generator) -> dict:180 M = int(self.providers.max()) + 1181 return {182 "c": rng.standard_normal(n_draws) * self.sigma_c,183 "u": rng.standard_normal((n_draws, self.J)) * self.sigma_u,184 "g": rng.standard_normal((n_draws, M)) * self.sigma_g,185 "e": rng.standard_normal((n_draws, self.K)),186 }187188189def exchange_views(views: list[np.ndarray], cited: list[np.ndarray] | None = None, J: int | None = None) -> list[np.ndarray]:190 """Evidence exchange: each agent adds the ids cited by the others (default: the union of all views)."""191 if cited is None:192 union = np.unique(np.concatenate([np.asarray(v) for v in views]))193 return [union.copy() for _ in views]194 return [np.unique(np.concatenate([np.asarray(v), np.asarray(c)])) for v, c in zip(views, cited)]195196197def random_partition_views(K: int, J: int, view_size: int, rng: np.random.Generator) -> list[np.ndarray]:198 """K views of `view_size` evidence ids each drawn without replacement (partially overlapping)."""199 return [np.sort(rng.choice(J, size=view_size, replace=False)) for _ in range(K)]200201202def placebo_mse(committee: EvidenceCommittee, perm: np.ndarray) -> float:203 """Equal-weight pooled MSE when the prompt-specific variances are permuted across slots (views fixed)."""204 return committee.mse_equal_closed_form(tau2=committee.tau2[perm])205206207def simulate_deliberation(committee: EvidenceCommittee, rounds: int, lam: float, n_draws: int,208 rng: np.random.Generator, anchor: str = "median", exchange_at: int | None = None,209 persuasion_sd: float = 0.0) -> dict[str, np.ndarray]:210 """Trajectory of dispersion and pooled errors under opinion-only deliberation, optional evidence211 exchange at round `exchange_at` (views enlarged to the union once) and optional persuasion noise212 (zero-mean revisions independent of the errors). Returns per-round arrays."""213 shocks = committee.draw_shocks(n_draws, rng)214 x = committee.draw(n_draws, rng, shocks=shocks)215 disp, mse_mean, mse_median = [], [], []216217 def record(x):218 disp.append(float(np.mean(x.std(axis=1, ddof=1))))219 mse_mean.append(float(np.mean(x.mean(axis=1) ** 2)))220 mse_median.append(float(np.mean(np.median(x, axis=1) ** 2)))221222 record(x)223 for r in range(1, rounds + 1):224 stat = np.median(x, axis=1) if anchor == "median" else x.mean(axis=1)225 x = (1.0 - lam) * x + lam * stat[:, None]226 if exchange_at is not None and r == exchange_at:227 # The exchange replaces each agent's evidence component by the enlarged-view component,228 # keeping the same evidence draws (common random numbers) and the same reasoning errors.229 A_old = committee.view_matrix()230 A_new = committee.view_matrix(exchange_views(committee.views))231 x = x + shocks["u"] @ (A_new - A_old).T232 if persuasion_sd > 0:233 x = x + rng.standard_normal(x.shape) * persuasion_sd234 record(x)235 return {"dispersion": np.array(disp), "mse_mean": np.array(mse_mean), "mse_median": np.array(mse_median)}236