//
//  Summarize.swift
//  Metrika
//
//  Author:  Simon-Pierre Boucher
//  Contact: contact@spboucher.ai
//  Copyright © 2026 Simon-Pierre Boucher. All rights reserved.
//

import Foundation

/// Descriptive statistics for one variable, Stata `summarize` semantics.
public struct ZQSummary: Equatable, Sendable {
    public var name: String
    public var observationCount: Int
    public var mean: Double
    /// Sample variance (n−1 denominator).
    public var variance: Double
    public var standardDeviation: Double
    public var minimum: Double
    public var maximum: Double
    /// Detail block (percentiles, skewness, kurtosis) — only populated for
    /// `summarize, detail`.
    public var detail: Detail?

    public struct Detail: Equatable, Sendable {
        /// Percentiles keyed by percent (1, 5, 10, 25, 50, 75, 90, 95, 99),
        /// computed with Stata's definition: for np/100 integral, the mean
        /// of adjacent order statistics; otherwise the next order statistic.
        public var percentiles: [Int: Double]
        /// Skewness m₃/m₂^{3/2} and kurtosis m₄/m₂² on central moments
        /// with n denominators (Stata definitions).
        public var skewness: Double
        public var kurtosis: Double
        public var smallest: [Double]
        public var largest: [Double]
    }
}

public enum ZQSummarize {

    /// Summarizes non-missing values. Returns nil when all values are
    /// missing (Stata prints a zero-observation row in that case; the
    /// caller renders it).
    public static func summary(
        name: String, values: [Double], missing: [Bool], detail: Bool = false
    ) -> ZQSummary? {
        var data: [Double] = []
        data.reserveCapacity(values.count)
        for i in 0..<values.count where !missing[i] {
            data.append(values[i])
        }
        let n = data.count
        guard n > 0 else { return nil }

        var sum = 0.0
        for v in data { sum += v }
        let mean = sum / Double(n)

        var m2 = 0.0, m3 = 0.0, m4 = 0.0
        var minimum = data[0], maximum = data[0]
        for v in data {
            let d = v - mean
            let d2 = d * d
            m2 += d2
            m3 += d2 * d
            m4 += d2 * d2
            if v < minimum { minimum = v }
            if v > maximum { maximum = v }
        }
        let variance = n > 1 ? m2 / Double(n - 1) : 0

        var summary = ZQSummary(
            name: name,
            observationCount: n,
            mean: mean,
            variance: variance,
            standardDeviation: variance.squareRoot(),
            minimum: minimum,
            maximum: maximum,
            detail: nil
        )

        if detail {
            let sorted = data.sorted()
            var percentiles: [Int: Double] = [:]
            for p in [1, 5, 10, 25, 50, 75, 90, 95, 99] {
                percentiles[p] = stataPercentile(sorted: sorted, percent: Double(p))
            }
            let mc2 = m2 / Double(n)
            let mc3 = m3 / Double(n)
            let mc4 = m4 / Double(n)
            summary.detail = ZQSummary.Detail(
                percentiles: percentiles,
                skewness: mc2 > 0 ? mc3 / pow(mc2, 1.5) : .nan,
                kurtosis: mc2 > 0 ? mc4 / (mc2 * mc2) : .nan,
                smallest: Array(sorted.prefix(4)),
                largest: Array(sorted.suffix(4))
            )
        }
        return summary
    }

    /// Stata percentile definition: with np = n·p/100, if np is integral
    /// the percentile is (x₍np₎ + x₍np+1₎)/2, otherwise x₍⌈np⌉₎ (1-based
    /// order statistics).
    public static func stataPercentile(sorted: [Double], percent: Double) -> Double {
        let n = sorted.count
        precondition(n > 0)
        let np = Double(n) * percent / 100
        let isIntegral = np == np.rounded()
        if isIntegral {
            let index = Int(np)
            if index >= n { return sorted[n - 1] }
            if index < 1 { return sorted[0] }
            return 0.5 * (sorted[index - 1] + sorted[index])
        }
        let index = min(n, max(1, Int(np.rounded(.up))))
        return sorted[index - 1]
    }
}