diff --git a/Libs/Estimator.csproj b/Libs/Estimator.csproj index afc2def..a85ba4c 100644 --- a/Libs/Estimator.csproj +++ b/Libs/Estimator.csproj @@ -6,7 +6,7 @@ True - 1.0.3 + 1.0.4 Statistics and performance metrics in trading, CAGR, Sharpe, MAE, MFE, and others. artemiusgreat indemos.com @@ -15,5 +15,9 @@ finance trading metrics forex performance-metrics stock-market stocks futures strategies backtestings MIT + + + + diff --git a/Libs/Estimators/CAGR/CAGR.cs b/Libs/Estimators/CAGR/CAGR.cs index ccade18..5c32e53 100644 --- a/Libs/Estimators/CAGR/CAGR.cs +++ b/Libs/Estimators/CAGR/CAGR.cs @@ -1,4 +1,3 @@ -using Estimator.Extensions; using Estimator.Models; using System; using System.Collections.Generic; @@ -21,6 +20,6 @@ public class CAGR /// /// Calculate /// - public virtual double Calculate() => (Math.Pow(Items.Sum(o => o.Value), 1.0 / Count) - 1.0).Round(); + public virtual double Calculate() => Math.Pow(Items.Sum(o => o.Value), 1.0 / Count) - 1.0; } } diff --git a/Libs/Estimators/Edge/EdgeRatio.cs b/Libs/Estimators/Edge/EdgeRatio.cs index 73b2285..c9c6136 100644 --- a/Libs/Estimators/Edge/EdgeRatio.cs +++ b/Libs/Estimators/Edge/EdgeRatio.cs @@ -1,5 +1,5 @@ -using Estimator.Extensions; using Estimator.Models; +using System; using System.Collections.Generic; namespace Estimator.Estimators @@ -19,7 +19,7 @@ public virtual double Calculate() var averageGain = new MFE { Items = Items }.Calculate(); var averageLoss = new MAE { Items = Items }.Calculate(); - return (averageGain / averageLoss).Round(); + return averageGain / averageLoss; } } } diff --git a/Libs/Estimators/Kestner/KestnerRatio.cs b/Libs/Estimators/Kestner/KestnerRatio.cs index df4e6fc..c156b85 100644 --- a/Libs/Estimators/Kestner/KestnerRatio.cs +++ b/Libs/Estimators/Kestner/KestnerRatio.cs @@ -1,4 +1,3 @@ -using Estimator.Extensions; using Estimator.Models; using System; using System.Collections.Generic; @@ -26,7 +25,7 @@ public virtual double Calculate() error = 1.0; } - return (slope / error).Round(); + return slope / error; } } } diff --git a/Libs/Estimators/MAE/MAE.cs b/Libs/Estimators/MAE/MAE.cs index aeff199..b6de6d7 100644 --- a/Libs/Estimators/MAE/MAE.cs +++ b/Libs/Estimators/MAE/MAE.cs @@ -1,5 +1,5 @@ -using Estimator.Extensions; using Estimator.Models; +using System; using System.Collections.Generic; using System.Linq; @@ -15,6 +15,6 @@ public class MAE /// /// Calculate /// - public virtual double Calculate() => Items.Average(o => o.Value - o.Min).Round(); + public virtual double Calculate() => Items.Average(o => o.Value - o.Min); } } diff --git a/Libs/Estimators/MFE/MFE.cs b/Libs/Estimators/MFE/MFE.cs index 2c6fed2..59f5ec7 100644 --- a/Libs/Estimators/MFE/MFE.cs +++ b/Libs/Estimators/MFE/MFE.cs @@ -1,5 +1,5 @@ -using Estimator.Extensions; using Estimator.Models; +using System; using System.Collections.Generic; using System.Linq; @@ -15,6 +15,6 @@ public class MFE /// /// Calculate /// - public virtual double Calculate() => Items.Average(o => o.Max - o.Value).Round(); + public virtual double Calculate() => Items.Average(o => o.Max - o.Value); } } diff --git a/Libs/Estimators/PF/PF.cs b/Libs/Estimators/PF/PF.cs index 3876169..e41b680 100644 --- a/Libs/Estimators/PF/PF.cs +++ b/Libs/Estimators/PF/PF.cs @@ -1,4 +1,3 @@ -using Estimator.Extensions; using Estimator.Models; using System; using System.Collections.Generic; @@ -26,7 +25,7 @@ public virtual double Calculate() gains = 1.0; } - return Math.Abs(gains / losses).Round(); + return Math.Abs(gains / losses); } } } diff --git a/Libs/Estimators/RAR/RAR.cs b/Libs/Estimators/RAR/RAR.cs index e7c41e1..2f3d816 100644 --- a/Libs/Estimators/RAR/RAR.cs +++ b/Libs/Estimators/RAR/RAR.cs @@ -1,4 +1,3 @@ -using Estimator.Extensions; using Estimator.Models; using System; using System.Collections.Generic; @@ -70,7 +69,7 @@ public virtual double Calculate() denominator = 1.0; } - return (Math.Sqrt(Count) * (mean - Rate) / denominator).Round(); + return Math.Sqrt(Count) * (mean - Rate) / denominator; } } } diff --git a/Libs/Estimators/Regression/Regression.cs b/Libs/Estimators/Regression/Regression.cs index a4c0766..25cd221 100644 --- a/Libs/Estimators/Regression/Regression.cs +++ b/Libs/Estimators/Regression/Regression.cs @@ -1,4 +1,3 @@ -using Estimator.Extensions; using Estimator.Models; using System; using System.Collections.Generic; diff --git a/Libs/Estimators/StandardScore/StandardScore.cs b/Libs/Estimators/StandardScore/StandardScore.cs index c4a101a..e09c98b 100644 --- a/Libs/Estimators/StandardScore/StandardScore.cs +++ b/Libs/Estimators/StandardScore/StandardScore.cs @@ -1,4 +1,3 @@ -using Estimator.Extensions; using Estimator.Models; using System; using System.Collections.Generic; @@ -48,7 +47,7 @@ public virtual double Calculate() deviation = 1.0; } - return ((seriesCount - mean - 0.5) / deviation).Round(); + return (seriesCount - mean - 0.5) / deviation; } } } diff --git a/Libs/Extensions/Double.cs b/Libs/Extensions/Double.cs deleted file mode 100644 index 1767b7c..0000000 --- a/Libs/Extensions/Double.cs +++ /dev/null @@ -1,12 +0,0 @@ -using System; - -namespace Estimator.Extensions -{ - public static class DoubleExtensions - { - public static double Round(this double input, int points = 2, MidpointRounding mode = MidpointRounding.ToEven) - { - return Math.Round(input, points, mode); - } - } -} diff --git a/Libs/Services/CointegrationService.cs b/Libs/Services/CointegrationService.cs new file mode 100644 index 0000000..d9d625d --- /dev/null +++ b/Libs/Services/CointegrationService.cs @@ -0,0 +1,147 @@ +using MathNet.Numerics; +using MathNet.Numerics.LinearAlgebra; +using MathNet.Numerics.LinearAlgebra.Factorization; +using System; +using System.Linq; + +namespace Estimator.Services +{ + /// + /// References + /// Numerically Stable Cointegration Analysis - Jurgen Doornik + /// Alternative asymptotics for cointegration tests in large VARs - Alexei Onatski and Chen Wang + /// Testing for Cointegration Using the Johansen Methodology when Variables are Near-Integrated - Erik Hjalmarsson + /// Time Series Analysis - Hamilton, J.D. + /// https://github.com/iisayoo/johansen/blob/master/johansen/johansen.py + /// + public class CointegrationService + { + /// + /// Cointegration + /// + /// + /// + /// + public static (double, Vector) Johansen(Matrix inputMatrix, int lags = 1) + { + var matrix = Center(inputMatrix); + var variables = matrix.ColumnCount; + var observations = matrix.RowCount - Math.Max(1, lags); + + // Step 1: First difference + var xSub = + matrix.SubMatrix(1, matrix.RowCount - 1, 0, variables) - + matrix.SubMatrix(0, matrix.RowCount - 1, 0, variables); + + // Step 2: Create lagged matrices + var xLag = Lag(matrix, lags); + var xSubLag = Lag(xSub, lags); + + // Step 3: Resize matrices + xSub = xSub.Resize(observations, xSub.ColumnCount); + xLag = xLag.Resize(observations, xLag.ColumnCount); + xSubLag = xSubLag.Resize(observations, xSubLag.ColumnCount); + + // Step 4: Compute residuals of xSub and xLag on xSubLag + var xSubLagInv = xSubLag.PseudoInverse(); + var u = xSub - xSubLag * xSubLagInv * xSub; + var v = xLag - xSubLag * xSubLagInv * xLag; + + // Step 5: Compute eigenvalues and eigenvectors + var Suu = u.TransposeThisAndMultiply(u) / observations; + var Svv = v.TransposeThisAndMultiply(v) / observations; + var Suv = u.TransposeThisAndMultiply(v) / observations; + var evd = (Svv.PseudoInverse() * Suv.TransposeThisAndMultiply(Suu.PseudoInverse() * Suv)).Evd(); + var eigenValues = evd.EigenValues.Real(); + var maxIndex = eigenValues.MaximumIndex(); + var eigenVectors = evd.EigenVectors.Column(maxIndex); + + // Step 6: Compare trace statistics with critical values to determine cointegration rank + var rank = GetRank(eigenValues, observations); + + return (rank, eigenVectors); + } + + /// + /// Subtract mean + /// + /// + /// + private static Matrix Center(Matrix matrix) + { + var columnMeans = matrix.ColumnSums() / matrix.RowCount; + var mean = Matrix.Build.Dense(matrix.RowCount, matrix.ColumnCount, (i, ii) => columnMeans[ii]); + return matrix - mean; + } + + /// + /// Get approximated critical value for specified rank + /// + /// + /// + /// + private static double GetCriticalValue(int rank, int numVariables) + { + // Higher base critical value for better scaling + double baseCriticalValue = 3.76; + + // Stronger scaling factors to match real-world growth in critical values + double variableFactor = Math.Pow(2.5, numVariables - 2); // More aggressive scaling for variables + double rankFactor = Math.Pow(2.5, rank); // Steep scaling for rank + + // Final critical value is the base value multiplied by both factors + double criticalValue = baseCriticalValue * variableFactor * rankFactor; + + return criticalValue; + } + + /// + /// Lag matrix + /// + /// + /// + /// + /// + private static Matrix Lag(Matrix matrix, int lags = 1) + { + var response = Matrix.Build.Dense(matrix.RowCount, matrix.ColumnCount * lags); + + for (var i = 0; i < lags; i++) + { + var subMatrix = matrix.SubMatrix(0, matrix.RowCount - i, 0, matrix.ColumnCount); + response.SetSubMatrix(i + 1, subMatrix.RowCount - 1, matrix.ColumnCount * i, subMatrix.ColumnCount, subMatrix); + } + + return response; + } + + /// + /// Get + /// + /// + /// + /// + private static int GetRank(Vector eigenValues, int observations) + { + var variables = eigenValues.Count; + + // Loop through from rank r to m (the number of variables) + for (var i = 0; i < variables; i++) + { + // Calculate the test statistic: -t * sum(log(1 - eigenvalues[i:])) + var statistic = -observations * eigenValues.SubVector(i + 1, variables - i - 1).Sum(o => Math.Log(1 - o)); + // Get the corresponding critical value for this rank + var critValue = GetCriticalValue(i, variables); + + // If the statistic is less than the critical value, return the rank i + if (statistic < critValue) + { + return i; + } + } + + // If no rank satisfies the condition, return m + return variables; + } + } +} diff --git a/Libs/Services/OptionService.cs b/Libs/Services/OptionService.cs new file mode 100644 index 0000000..7775e78 --- /dev/null +++ b/Libs/Services/OptionService.cs @@ -0,0 +1,257 @@ +using MathNet.Numerics.Distributions; +using MathNet.Numerics.RootFinding; +using System; + +namespace Estimator.Services +{ + public enum OptionSideEnum : byte + { + None = 0, + Put = 1, + Call = 2, + Share = 3 + } + + public class OptionService + { + /// + /// Asset or nothing + /// + /// + /// + /// + /// + /// + /// + /// + private static double D1(double S, double K, double T, double sigma, double r, double q) + { + return (Math.Log(S / K) + (r - q + (sigma * sigma) / 2) * T) / (sigma * Math.Sqrt(T)); + } + + /// + /// Cash or nothing + /// + /// + /// + /// + /// + private static double D2(double T, double sigma, double d1) + { + return d1 - sigma * Math.Sqrt(T); + } + + /// + /// Computes theoretical price. + /// + /// call or put + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Premium(OptionSideEnum optionType, double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + double d2 = D2(T, sigma, d1); + + switch (optionType) + { + case OptionSideEnum.Call: + return S * Math.Exp(-q * T) * Normal.CDF(0, 1, d1) - K * Math.Exp(-r * T) * Normal.CDF(0, 1, d2); + + case OptionSideEnum.Put: + return K * Math.Exp(-r * T) * Normal.CDF(0, 1, -d2) - S * Math.Exp(-q * T) * Normal.CDF(0, 1, -d1); + } + + return 0; + } + + /// + /// Computes Vega. The amount of option price change for each 1% change in vol (sigma) + /// + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Vega(double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + double vega = S * Math.Exp(-q * T) * Normal.PDF(0, 1, d1) * Math.Sqrt(T); + return vega / 100; + } + + /// + /// Volatility + /// + /// + /// + /// + /// + /// + /// + /// + /// + public static double IV(OptionSideEnum optionType, double S, double K, double T, double r, double q, double optionMarketPrice) + { + Func f = sigma => Premium(optionType, S, K, T, sigma, r, q) - optionMarketPrice; + Func df = sigma => Vega(S, K, T, sigma, r, q); + + return RobustNewtonRaphson.FindRoot(f, df, lowerBound: 0, upperBound: 100, accuracy: 0.001); + } + + /// + /// Computes theta. + /// + /// call or put + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Theta(OptionSideEnum optionType, double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + double d2 = D2(T, sigma, d1); + + switch (optionType) + { + case OptionSideEnum.Call: + { + double theta = -Math.Exp(-q * T) * (S * Normal.PDF(0, 1, d1) * sigma) / (2.0 * Math.Sqrt(T)) + - (r * K * Math.Exp(-r * T) * Normal.CDF(0, 1, d2)) + + q * S * Math.Exp(-q * T) * Normal.CDF(0, 1, d1); + + return theta / 365; + } + + case OptionSideEnum.Put: + { + double theta = -Math.Exp(-q * T) * (S * Normal.PDF(0, 1, d1) * sigma) / (2.0 * Math.Sqrt(T)) + + (r * K * Math.Exp(-r * T) * Normal.PDF(0, 1, -d2)) + - q * S * Math.Exp(-q * T) * Normal.CDF(0, 1, -d1); + + return theta / 365; + } + + default: + throw new NotSupportedException(); + } + } + + /// + /// Computes delta. + /// + /// call or put + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Delta(OptionSideEnum optionType, double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + + switch (optionType) + { + case OptionSideEnum.Call: + return Math.Exp(-r * T) * Normal.CDF(0, 1, d1); + + case OptionSideEnum.Put: + return -Math.Exp(-r * T) * Normal.CDF(0, 1, -d1); + + default: + throw new NotSupportedException(); + } + } + + /// + /// Computes gamma. + /// + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Gamma(double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + return Math.Exp(-q * T) * (Normal.PDF(0, 1, d1) / (S * sigma * Math.Sqrt(T))); + } + + /// + /// Computes delta. + /// + /// call or put + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Rho(OptionSideEnum optionType, double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + double d2 = D2(T, sigma, d1); + + switch (optionType) + { + case OptionSideEnum.Call: + return K * T * Math.Exp(-r * T) * Normal.CDF(0, 1, d2); + + case OptionSideEnum.Put: + return -K * T * Math.Exp(-r * T) * Normal.CDF(0, 1, -d2); + + default: + throw new NotSupportedException(); + } + } + + /// + /// Computes Vanna. The amount of option price change for each 1% change in vol (sigma) + /// + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Vanna(double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + double vega = Vega(S, K, T, sigma, r, q); + return (d1 / sigma) * vega; + } + + /// + /// Computes Volga. The amount of option price change for each 1% change in vol (sigma) + /// + /// Underlying price + /// Strike price + /// Time to expiration in % of year + /// Volatility + /// continuously compounded risk-free interest rate + /// continuously compounded dividend yield + /// + public static double Volga(double S, double K, double T, double sigma, double r, double q) + { + double d1 = D1(S, K, T, sigma, r, q); + double d2 = D2(T, sigma, d1); + double vega = Vega(S, K, T, sigma, r, q); + return vega * d1 * d2 / sigma; + } + } +}