Added tail dependency to copulas. Added seasonal decomposition to times series class. Fixed some errors in docs.

Author: HadenSmith

Numerics/Data/Time Series/TimeSeries.cs

@@ -30,6 +30,7 @@
 
 using Numerics.Data.Statistics;
 using Numerics.MachineLearning;
+using Numerics.Mathematics;
 using Numerics.Sampling;
 using System;
 using System.Collections.Generic;
@@ -1019,6 +1020,101 @@ public TimeSeries MovingSum(int period)
             return timeSeries;
         }
 
+        /// <summary>
+        /// Performs classical additive seasonal decomposition using FFT-based seasonal extraction.
+        /// </summary>
+        /// <param name="period">The seasonal period (e.g., 12 for monthly data with annual seasonality).</param>
+        /// <returns>
+        /// A tuple containing:
+        /// <list type="bullet">
+        /// <item><description>Trend: The trend component as a TimeSeries (moving average with the given period).</description></item>
+        /// <item><description>Seasonal: The seasonal component as a double array of length equal to the original series.</description></item>
+        /// <item><description>Residual: The residual component as a TimeSeries (defined where the trend is defined).</description></item>
+        /// </list>
+        /// </returns>
+        /// <exception cref="ArgumentException">Thrown when the period is less than 2 or the series contains fewer than 2 complete periods.</exception>
+        public (TimeSeries Trend, double[] Seasonal, TimeSeries Residual) SeasonalDecompose(int period)
+        {
+            if (period < 2)
+                throw new ArgumentException("Period must be at least 2.", nameof(period));
+            if (Count < 2 * period)
+                throw new ArgumentException("Time series must contain at least 2 complete periods.", nameof(period));
+
+            SortByTime();
+            int n = Count;
+
+            // Step 1: Compute trend via moving average
+            var trend = MovingAverage(period);
+
+            // Step 2: Build full-length arrays and detrend
+            double[] values = new double[n];
+            double[] trendFull = new double[n];
+            bool[] hasTrend = new bool[n];
+
+            for (int i = 0; i < n; i++)
+                values[i] = this[i].Value;
+
+            // Map trend values to original indices
+            // MovingAverage outputs (Count - period + 1) values starting at index (period - 1)
+            int trendStart = period - 1;
+            for (int i = 0; i < trend.Count; i++)
+            {
+                trendFull[trendStart + i] = trend[i].Value;
+                hasTrend[trendStart + i] = true;
+            }
+
+            // Create detrended series
+            double[] detrended = new double[n];
+            for (int i = 0; i < n; i++)
+                detrended[i] = hasTrend[i] ? values[i] - trendFull[i] : 0.0;
+
+            // Step 3: FFT-based seasonal extraction
+            // Pad to power of 2 for FFT
+            int fftLength = (int)Math.Pow(2, Math.Ceiling(Math.Log(n, 2)));
+            if (fftLength < n) fftLength *= 2;
+            double[] fftData = new double[fftLength];
+            Array.Copy(detrended, fftData, n);
+
+            // Forward FFT
+            Fourier.RealFFT(fftData);
+
+            // Identify harmonic bins of the seasonal frequency
+            // Harmonic h corresponds to frequency bin k = round(h * fftLength / period)
+            var harmonicBins = new HashSet<int>();
+            for (int h = 1; ; h++)
+            {
+                int k = (int)Math.Round((double)h * fftLength / period);
+                if (k <= 0 || k >= fftLength / 2) break;
+                harmonicBins.Add(k);
+            }
+
+            // Keep only harmonic frequencies
+            double[] filtered = new double[fftLength];
+            foreach (int k in harmonicBins)
+            {
+                filtered[2 * k] = fftData[2 * k];
+                filtered[2 * k + 1] = fftData[2 * k + 1];
+            }
+
+            // Inverse FFT
+            Fourier.RealFFT(filtered, true);
+
+            // Scale by 2/fftLength as required by the RealFFT inverse transform
+            double[] seasonal = new double[n];
+            for (int i = 0; i < n; i++)
+                seasonal[i] = filtered[i] * 2.0 / fftLength;
+
+            // Step 4: Compute residual
+            var residual = new TimeSeries(TimeInterval);
+            for (int i = 0; i < n; i++)
+            {
+                if (hasTrend[i])
+                    residual.Add(new SeriesOrdinate<DateTime, double>(this[i].Index, values[i] - trendFull[i] - seasonal[i]));
+            }
+
+            return (trend, seasonal, residual);
+        }
+
         /// <summary>
         /// Shift all of the dates to match the new start date.
         /// </summary>

Numerics/Distributions/Bivariate Copulas/AMHCopula.cs

@@ -182,6 +182,18 @@ public override double[] InverseCDF(double u, double v)
             return [u, v];
         }
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U = 0.
+        /// The AMH copula has no tail dependence.
+        /// </summary>
+        public override double UpperTailDependence => 0.0;
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L = 0.
+        /// The AMH copula has no tail dependence.
+        /// </summary>
+        public override double LowerTailDependence => 0.0;
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {

Numerics/Distributions/Bivariate Copulas/Base/BivariateCopula.cs

@@ -119,6 +119,12 @@ public bool ParametersValid
         /// <inheritdoc/>
         public abstract string ShortDisplayName { get; }
 
+        /// <inheritdoc/>
+        public abstract double UpperTailDependence { get; }
+
+        /// <inheritdoc/>
+        public abstract double LowerTailDependence { get; }
+
         #endregion
 
         #region Methods

Numerics/Distributions/Bivariate Copulas/Base/IBivariateCopula.cs

@@ -142,6 +142,16 @@ public interface IBivariateCopula : IDistribution
         /// <param name="v">Probability between 0 and 1.</param>
         double[] InverseCDF(double u, double v);
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U.
+        /// </summary>
+        double UpperTailDependence { get; }
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L.
+        /// </summary>
+        double LowerTailDependence { get; }
+
         /// <summary>
         /// Generate random values of a distribution given a sample size.
         /// </summary>

Numerics/Distributions/Bivariate Copulas/ClaytonCopula.cs

@@ -160,6 +160,24 @@ public override double[] InverseCDF(double u, double v)
             return [u, v];
         }
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U = 0.
+        /// The Clayton copula has no upper tail dependence.
+        /// </summary>
+        public override double UpperTailDependence => 0.0;
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L = 2^(-1/θ).
+        /// </summary>
+        public override double LowerTailDependence
+        {
+            get
+            {
+                if (Theta <= 0.0) return 0.0;
+                return Math.Pow(2.0, -1.0 / Theta);
+            }
+        }
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {

Numerics/Distributions/Bivariate Copulas/FrankCopula.cs

@@ -183,6 +183,18 @@ public override double[] InverseCDF(double u, double v)
             return [u, v];
         }
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U = 0.
+        /// The Frank copula has no tail dependence.
+        /// </summary>
+        public override double UpperTailDependence => 0.0;
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L = 0.
+        /// The Frank copula has no tail dependence.
+        /// </summary>
+        public override double LowerTailDependence => 0.0;
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {

Numerics/Distributions/Bivariate Copulas/GumbelCopula.cs

@@ -159,6 +159,23 @@ public override double[] InverseCDF(double u, double v)
             return [u, v];
         }
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U = 2 - 2^(1/θ).
+        /// </summary>
+        public override double UpperTailDependence
+        {
+            get
+            {
+                return 2.0 - Math.Pow(2.0, 1.0 / Theta);
+            }
+        }
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L = 0.
+        /// The Gumbel copula has no lower tail dependence.
+        /// </summary>
+        public override double LowerTailDependence => 0.0;
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {

Numerics/Distributions/Bivariate Copulas/JoeCopula.cs

@@ -161,6 +161,23 @@ public override double[] InverseCDF(double u, double v)
             return [u, v];
         }
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U = 2 - 2^(1/θ).
+        /// </summary>
+        public override double UpperTailDependence
+        {
+            get
+            {
+                return 2.0 - Math.Pow(2.0, 1.0 / Theta);
+            }
+        }
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L = 0.
+        /// The Joe copula has no lower tail dependence.
+        /// </summary>
+        public override double LowerTailDependence => 0.0;
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {

Numerics/Distributions/Bivariate Copulas/NormalCopula.cs

@@ -190,6 +190,18 @@ public override double[] InverseCDF(double u, double v)
             return [u, v];
         }
 
+        /// <summary>
+        /// Gets the upper tail dependence coefficient λ_U = 0.
+        /// The Normal copula has no upper tail dependence.
+        /// </summary>
+        public override double UpperTailDependence => 0.0;
+
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L = 0.
+        /// The Normal copula has no lower tail dependence.
+        /// </summary>
+        public override double LowerTailDependence => 0.0;
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {

Numerics/Distributions/Bivariate Copulas/StudentTCopula.cs

@@ -365,19 +365,19 @@ public override double[] InverseCDF(double u, double v)
         }
 
         /// <summary>
-        /// Gets the symmetric tail dependence coefficient λ = λ_L = λ_U.
+        /// Gets the upper tail dependence coefficient λ_U.
         /// </summary>
         /// <remarks>
         /// <para>
         /// The t-copula has symmetric upper and lower tail dependence:
         /// <code>
-        ///     λ = 2 · t_{ν+1}(-√((ν+1)(1-ρ)/(1+ρ)))
+        ///     λ_U = λ_L = 2 · t_{ν+1}(-√((ν+1)(1-ρ)/(1+ρ)))
         /// </code>
         /// where t_{ν+1} is the CDF of the univariate Student's t with ν+1 degrees of freedom.
         /// For ρ = -1, λ = 0. For ρ = 1, λ = 1. As ν → ∞, λ → 0 (Normal copula limit).
         /// </para>
         /// </remarks>
-        public double TailDependence
+        public override double UpperTailDependence
         {
             get
             {
@@ -393,6 +393,12 @@ public double TailDependence
             }
         }
 
+        /// <summary>
+        /// Gets the lower tail dependence coefficient λ_L.
+        /// The t-copula has symmetric tail dependence, so λ_L = λ_U.
+        /// </summary>
+        public override double LowerTailDependence => UpperTailDependence;
+
         /// <inheritdoc/>
         public override BivariateCopula Clone()
         {