feat(football): implement World Cup 2026 simulator and Empirical-Poisson Mixture Model

pkg/football/logic/score_model.go

@@ -0,0 +1,191 @@
+package logic
+
+import (
+	"math"
+	"math/rand"
+	"sync"
+)
+
+// RawScoreRecord represents the raw statistics for a scoreline in football history.
+type RawScoreRecord struct {
+	GoalsA int `json:"goals_a"`
+	GoalsB int `json:"goals_b"`
+	Count  int `json:"count"`
+}
+
+// WorldCupRawScores represents the actual historical frequencies of scorelines
+// in the FIFA World Cup history (964 matches from 1930 to 2022).
+var WorldCupRawScores = []RawScoreRecord{
+	{GoalsA: 1, GoalsB: 0, Count: 182},
+	{GoalsA: 2, GoalsB: 1, Count: 152},
+	{GoalsA: 2, GoalsB: 0, Count: 111},
+	{GoalsA: 1, GoalsB: 1, Count: 92},
+	{GoalsA: 0, GoalsB: 0, Count: 78},
+	{GoalsA: 3, GoalsB: 1, Count: 68},
+	{GoalsA: 3, GoalsB: 0, Count: 57},
+	{GoalsA: 3, GoalsB: 2, Count: 43},
+	{GoalsA: 2, GoalsB: 2, Count: 35},
+	{GoalsA: 4, GoalsB: 1, Count: 31},
+	{GoalsA: 4, GoalsB: 0, Count: 24},
+	{GoalsA: 4, GoalsB: 2, Count: 17},
+	{GoalsA: 6, GoalsB: 1, Count: 11},
+	{GoalsA: 5, GoalsB: 2, Count: 9},
+	{GoalsA: 5, GoalsB: 0, Count: 7},
+	{GoalsA: 3, GoalsB: 3, Count: 7},
+	{GoalsA: 5, GoalsB: 1, Count: 7},
+	{GoalsA: 6, GoalsB: 0, Count: 5},
+	{GoalsA: 7, GoalsB: 0, Count: 5},
+	{GoalsA: 4, GoalsB: 3, Count: 3},
+	{GoalsA: 7, GoalsB: 1, Count: 3},
+	{GoalsA: 8, GoalsB: 1, Count: 3},
+	{GoalsA: 4, GoalsB: 4, Count: 2},
+	{GoalsA: 6, GoalsB: 3, Count: 2},
+	{GoalsA: 9, GoalsB: 0, Count: 2},
+	{GoalsA: 5, GoalsB: 3, Count: 1},
+	{GoalsA: 6, GoalsB: 2, Count: 1},
+	{GoalsA: 7, GoalsB: 2, Count: 1},
+	{GoalsA: 7, GoalsB: 3, Count: 1},
+	{GoalsA: 6, GoalsB: 5, Count: 1},
+	{GoalsA: 8, GoalsB: 3, Count: 1},
+	{GoalsA: 10, GoalsB: 1, Count: 1},
+	{GoalsA: 7, GoalsB: 5, Count: 1},
+}
+
+// poissonProbability computes the Poisson probability P(k; lambda).
+func poissonProbability(lambda float64, k int) float64 {
+	factorial := 1.0
+	for i := 1; i <= k; i++ {
+		factorial *= float64(i)
+	}
+	return (math.Pow(lambda, float64(k)) * math.Exp(-lambda)) / factorial
+}
+
+// DefaultScoreOutcomes returns the symmetric default probability distribution based on WorldCupRawScores.
+func DefaultScoreOutcomes() []ScoreOutcome {
+	var totalGoals int
+	var totalMatches int
+	for _, r := range WorldCupRawScores {
+		totalGoals += (r.GoalsA + r.GoalsB) * r.Count
+		totalMatches += r.Count
+	}
+
+	if totalMatches == 0 {
+		return []ScoreOutcome{{ScoreA: 1, ScoreB: 1, Probability: 1.0}}
+	}
+
+	// 1. Calculate expected goals per team (lambda)
+	lambda := float64(totalGoals) / (2.0 * float64(totalMatches))
+
+	// 2. Compute Poisson 1D and 2D probabilities (up to MaxGoals = 10)
+	maxGoals := 10
+	poisson1D := make([]float64, maxGoals+1)
+	var sumP float64
+	for g := 0; g <= maxGoals; g++ {
+		p := poissonProbability(lambda, g)
+		poisson1D[g] = p
+		sumP += p
+	}
+	for g := 0; g <= maxGoals; g++ {
+		poisson1D[g] /= sumP
+	}
+
+	// 3. Build Empirical 2D map
+	type scoreKey struct {
+		goalsA int
+		goalsB int
+	}
+	empirical := make(map[scoreKey]float64)
+	for _, r := range WorldCupRawScores {
+		if r.GoalsA <= maxGoals && r.GoalsB <= maxGoals {
+			prob := float64(r.Count) / float64(totalMatches)
+			if r.GoalsA == r.GoalsB {
+				empirical[scoreKey{r.GoalsA, r.GoalsB}] = prob
+			} else {
+				empirical[scoreKey{r.GoalsA, r.GoalsB}] = prob / 2.0
+				empirical[scoreKey{r.GoalsB, r.GoalsA}] = prob / 2.0
+			}
+		}
+	}
+
+	// 4. Combine into Mixture Model (alpha = 0.99)
+	alpha := 0.99
+	var outcomes []ScoreOutcome
+	for ga := 0; ga <= maxGoals; ga++ {
+		for gb := 0; gb <= maxGoals; gb++ {
+			pPoisson := poisson1D[ga] * poisson1D[gb]
+			pEmpirical := empirical[scoreKey{ga, gb}]
+			pMix := alpha*pEmpirical + (1.0-alpha)*pPoisson
+
+			outcomes = append(outcomes, ScoreOutcome{
+				ScoreA:      ga,
+				ScoreB:      gb,
+				Probability: pMix,
+			})
+		}
+	}
+
+	// Normalize just in case of float64 precision drift
+	var sumMix float64
+	for _, o := range outcomes {
+		sumMix += o.Probability
+	}
+	for i := range outcomes {
+		outcomes[i].Probability /= sumMix
+	}
+
+	return outcomes
+}
+
+// ScoreModel holds the probability distribution used to sample random scores.
+type ScoreModel struct {
+	Outcomes []ScoreOutcome
+	mu       sync.RWMutex
+}
+
+// NewScoreModel creates a ScoreModel from custom score outcomes.
+// If customOutcomes is empty, it uses DefaultScoreOutcomes().
+func NewScoreModel(customOutcomes []ScoreOutcome) *ScoreModel {
+	sm := &ScoreModel{}
+	if len(customOutcomes) > 0 {
+		// Normalize custom outcomes to ensure they sum to 1.0
+		var sum float64
+		for _, o := range customOutcomes {
+			sum += o.Probability
+		}
+		if sum > 0 {
+			normalized := make([]ScoreOutcome, len(customOutcomes))
+			for i, o := range customOutcomes {
+				normalized[i] = ScoreOutcome{
+					ScoreA:      o.ScoreA,
+					ScoreB:      o.ScoreB,
+					Probability: o.Probability / sum,
+				}
+			}
+			sm.Outcomes = normalized
+			return sm
+		}
+	}
+	sm.Outcomes = DefaultScoreOutcomes()
+	return sm
+}
+
+// RandomScore samples a random score outcome based on the probability distribution.
+func (sm *ScoreModel) RandomScore(r *rand.Rand) ScoreOutcome {
+	sm.mu.RLock()
+	outcomes := sm.Outcomes
+	sm.mu.RUnlock()
+
+	if len(outcomes) == 0 {
+		return ScoreOutcome{ScoreA: 0, ScoreB: 0, Probability: 1.0}
+	}
+
+	p := r.Float64()
+	var cumulative float64
+	for _, outcome := range outcomes {
+		cumulative += outcome.Probability
+		if p <= cumulative {
+			return outcome
+		}
+	}
+	return outcomes[len(outcomes)-1]
+}