**Approach (brief)** I sample weight vectors and per-option/criterion scores, compute weighted totals, and tally wins. For weights I support Dirichlet (for correlated simplex) or independent per-criterion samplers. For scores I support parametric samplers (normal, beta) per option/criterion. I also record samples where the top-ranked option changes to inspect flips.**Code (example implementation)**python
import numpy as np
from collections import defaultdict, Counter
def draw_weights(weight_distributions):
# weight_distributions: dict criterion -> {"dirichlet": alpha} or {"dist": callable}
if any('dirichlet' in v for v in weight_distributions.values()):
# support single global dirichlet via provided alphas under key '__dirichlet__'
alpha = weight_distributions.get('__dirichlet__')
if alpha is not None:
return np.random.dirichlet(alpha)
# independent weights, normalize
samples = np.array([v['dist']() for v in weight_distributions.values()])
return samples / samples.sum()
def draw_scores(score_distributions, options, criteria):
s = {}
for o in options:
s[o] = {}
for c in criteria:
spec = score_distributions.get(o, {}).get(c)
if spec is None: # fallback to nominal
s[o][c] = None
else:
if spec['type']=='normal':
s[o][c] = np.random.normal(spec['mu'], spec['sigma'])
elif spec['type']=='beta':
s[o][c] = np.random.beta(spec['a'], spec['b'])
else:
s[o][c] = spec['sample']()
return s
def run_monte_carlo(scoring_matrix, weight_distributions, score_distributions, samples=10000):
options = list(scoring_matrix.keys())
criteria = list(next(iter(scoring_matrix.values())).keys())
wins = Counter()
total_scores = {o: [] for o in options}
flip_examples = []
for _ in range(samples):
w = draw_weights(weight_distributions)
s = draw_scores(score_distributions, options, criteria)
totals = {}
for i,o in enumerate(options):
tot = 0.0
for j,c in enumerate(criteria):
val = s[o][c] if s[o][c] is not None else scoring_matrix[o][c]
tot += w[j]*val
totals[o]=tot
total_scores[o].append(tot)
winner = max(totals, key=totals.get)
wins[winner]+=1
# record flips: if nominal winner different than sampled winner
nominal_totals = {o: sum(scoring_matrix[o][c] for c in criteria)}
# simple check: note one example if differs
if winner != max(nominal_totals, key=nominal_totals.get) and len(flip_examples)<50:
flip_examples.append({'weights':w.copy(),'scores':s,'winner':winner})
results = {
'win_prob': {o: wins[o]/samples for o in options},
'score_ci': {o: (np.percentile(total_scores[o],2.5), np.percentile(total_scores[o],97.5)) for o in options},
'flip_examples': flip_examples
}
return results
