C# Math Library
Updated: April 14th, 2015 (07:11 AM) Views / Replies: 15,938 / 15
Created: March 18th, 2010 (10:30 AM) by MXASJ Attachments: 1

(If you already have an account, login at the top of the page)
futures io is the largest futures trading community on the planet, with over 90,000 members. At

futures io , our goal has always been and always will be to create a friendly, positive, forward-thinking community where members can openly share and discuss everything the world of trading has to offer. The community is one of the friendliest you will find on any subject, with members going out of their way to help others. Some of the primary differences between

futures io and other trading sites revolve around the standards of our community. Those standards include a code of conduct for our members, as well as extremely high standards that govern which partners we do business with, and which products or services we recommend to our members.

At

futures io , our focus is on quality education. No hype, gimmicks, or secret sauce. The truth is: trading is hard. To succeed, you need to surround yourself with the right support system, educational content, and trading mentors – all of which you can find on

futures io , utilizing our social trading environment.

With

futures io , you can find

honest trading reviews on brokers, trading rooms, indicator packages, trading strategies, and much more . Our trading review process is highly moderated to ensure that only genuine users are allowed, so you don’t need to worry about fake reviews.

We are fundamentally different than most other trading sites:

We are here to help. Just let us know what you need.
We work extremely hard to keep things positive in our community.
We do not tolerate rude behavior, trolling, or vendors advertising in posts.
We firmly believe in and encourage sharing. The holy grail is within you, we can help you find it.
We expect our members to participate and become a part of the community. Help yourself by helping others.
You'll need to

register in order to view the content of the threads and start contributing to our community.

It's free and simple.
-- Big Mike, Site Administrator

January 9th, 2012, 08:22 AM
#11 (permalink )
Elite Member

Seattle, WA

Futures Experience: Intermediate

Platform: NinjaTrader

Broker/Data: Interactive Brokers/IQFeed

Favorite Futures: Stocks

Posts: 27 since Nov 2009

Thanks: 10 given,
66
received

Here's the ALGLIB file. It's a NinjaTrader "Indicator", but it's just the entire ALGLIB 3.4 C# math library merged into a single file, put into a namespace called "ALGLIB". It's callable from strategies also.
I don't have time right now to give much of an example, but below are some static examples of matrix operations (least-squares fit, matrix multiplication, inverse, pseudoinverse, left divide & right divide, an OLS routine originally from MATLAB, and some helper routines). This (static) approach is easy to use (no instantiation required); it's fine for indicators but it falls down if you use it in a strategy which you run through the optimizer. That's because "new" is computationally expensive in C#, if you do "new" repeatedly it hits your performance and fills up all your memory. The alternative is to instantiate classes which allocate storage once and then re-use it, which helps a lot but it's more complicated and you still can eventually hit a wall for other reasons.
The code below is a bunch of copy & pastes; I think everything needed to compile it is here, but I haven't tested to be sure.
To use it:
Code
using ALGLIB;
public static void lsfitlinear(double[] y, double[,] x, out double[] beta)
{
int nobs = y.Length;
if (nobs != x.GetLength(0))
throw new ArgumentException("ols: y and x must have same length");
int nvar = x.GetLength(1);
int info;
alglib.lsfitreport rep;
alglib.lsfitlinear(y, x, nobs, nvar, out info, out beta, out rep);
}
public static double[,] mtimes(double[,] a, double[,] b, bool transposeA, bool transposeB)
{
int ra = a.GetLength(transposeA ? 1 : 0);
int ca = a.GetLength(transposeA ? 0 : 1);
int cb = b.GetLength(transposeB ? 0 : 1);
if (ca != b.GetLength(transposeB ? 1 : 0))
throw new ArgumentException("mtimes: columns(matrix a) must equal rows(matrix b)");
double [,] output = new double[ra,cb];
alglib.ablas.rmatrixgemm(ra,cb,ca,1,a,0,0,transposeA ? 1 : 0,b,0,0,transposeB ? 1 : 0,0,ref output,0,0);
return output;
}
public static double[,] inverse(double[,] a)
{
if (a.GetLength(0) != a.GetLength(1))
throw new ArgumentException("inverse: matrix must be square");
double[,] output = copy(a);
alglib.matinv.matinvreport rep = new alglib.matinv.matinvreport();
int info = new int();
alglib.matinv.rmatrixinverse(ref output,output.GetLength(0),ref info,rep);
return output;
}
public static double[,] copy(double[,] a)
{ // create a copy of a matrix
int r = a.GetLength(0);
int c = a.GetLength(1);
double[,] output = new double[r,c];
for (int i=0; i < r; i++)
{
for (int j=0; j < c; j++)
{
output[i,j] = a[i,j];
}
}
return output;
}
public static double[,] pseudoinverse(double[,] a)
{
int r = a.GetLength(0);
int c = a.GetLength(1);
int minrc = Math.Min(r,c);
double[] w = new double[minrc];
double[,] winv = new double[minrc,minrc];
double[,] u = new double[r,minrc];
double[,] vt = new double[minrc,c];
double[,] temp = new double[c,minrc];
double[,] output = new double[c,r];
alglib.svd.rmatrixsvd(a,r,c,1,1,2,ref w,ref u, ref vt);
// generate inverse winv of matrix W from the vector w from SVD
// to do this we transpose the matrix, and take the inverse of each vector element
// see: http://en.wikipedia.org/wiki/Singular_value_decomposition
// http://www.alglib.net/matrixops/general/svd.php
for (int i = 0; i < minrc; i++)
{
for (int j = 0; j < minrc; j++)
{
winv[i,j] = ((i == j) && (w[i] > alglib.math.machineepsilon)) ? 1.0/w[i] : 0.0; // inverse of w
}
}
alglib.ablas.rmatrixgemm(c,minrc,minrc,1,vt,0,0,1,winv,0,0,0,0,ref temp,0,0);
alglib.ablas.rmatrixgemm(c,r,minrc,1,temp,0,0,0,u,0,0,1,0,ref output,0,0);
return output;
}
public static double[,] mldivide(double[,] a, double[,] b)
{ // matrix left divide: inv(a)*b
int r = a.GetLength(0);
if (r != b.GetLength(0))
throw new ArgumentException ("mldivide: both matrices must have the same number of rows");
int ca = a.GetLength(1);
int cb = b.GetLength(1);
double[,] output = new double[ca,cb];
alglib.ablas.rmatrixgemm(ca,cb,r,1,pseudoinverse(a),0,0,0,b,0,0,0,0,ref output,0,0);
return output;
}
public static double[,] mrdivide(double[,] b, double[,] a)
{ // matrix right divide: b*inv(a)
int c = a.GetLength(1);
if (c != b.GetLength(1))
throw new ArgumentException ("mrdivide: both matrices must have the same number of columns");
int ra = a.GetLength(0);
int rb = b.GetLength(0);
double[,] output = new double[rb,ra];
alglib.ablas.rmatrixgemm(rb,ra,c,1,b,0,0,0,pseudoinverse(a),0,0,0,0,ref output,0,0);
return output;
}
public static void ols(double[] y, double[,] x, out double[] beta, out double[] yhat, out double[] residuals, out double sige, out double rsqr, out double rbar /*, out double[] tstats, out double[] bstd, out double dw*/)
{
// function results=ols(y,x)
// % PURPOSE: least-squares regression
// %---------------------------------------------------
// % USAGE: results = ols(y,x)
// % where: y = dependent variable vector (nobs x 1)
// % x = independent variables matrix (nobs x nvar)
// %---------------------------------------------------
// % RETURNS: a structure
// % results.meth = 'ols'
// % results.beta = bhat (nvar x 1)
// % results.tstat = t-stats (nvar x 1)
// % results.bstd = std deviations for bhat (nvar x 1)
// % results.yhat = yhat (nobs x 1)
// % results.resid = residuals (nobs x 1)
// % results.sige = e'*e/(n-k) scalar
// % results.rsqr = rsquared scalar
// % results.rbar = rbar-squared scalar
// % results.dw = Durbin-Watson Statistic
// % results.nobs = nobs
// % results.nvar = nvars
// % results.y = y data vector (nobs x 1)
// % results.bint = (nvar x2 ) vector with 95% confidence intervals on beta
// %---------------------------------------------------
int nobs = y.Length;
if (nobs != x.GetLength(0))
throw new ArgumentException("ols: y and x must have same length");
int nvar = x.GetLength(1);
int info;
yhat = new double[nobs];
alglib.lsfitreport rep;
alglib.lsfitlinear(y,x,nobs,nvar,out info, out beta, out rep); // least squares fit
alglib.ablas.rmatrixmv(x.GetLength(0),x.GetLength(1),x,0,0,0,beta,0,ref yhat,0); // calculate yhat (predicted values of Y)
residuals = minus(y,yhat);
double sigu = mtimes(residuals,residuals);
sige = sigu/System.Convert.ToDouble(nobs-nvar);
double[] ym = minus(y,mean(y));
double rsqr1 = sigu;
double rsqr2 = mtimes(ym,ym);
rsqr = 1.0 - rsqr1/rsqr2; // r-squared
rsqr1 = rsqr1/System.Convert.ToDouble(nobs-nvar);
rsqr2 = rsqr2/System.Convert.ToDouble(nobs-1);
rbar = rsqr2 == 0 ? rsqr : 1.0 - rsqr1/rsqr2;
// bstd, dw, tstat, bint not implemented
}
public static double mtimes(double[] x, double[] y)
{ // multiply 2 vectors (matrix multiply on 1-row and 1-col matrices)
int lx = x.Length;
if (lx != y.Length)
throw new ArgumentException("mtimes - input vectors must have the same length");
double output = 0;
for (int i = 0; i < lx; i++)
output += x[i] * y[i];
return output;
}
public static double mean(double[] x)
{ // mean of a vector
int len = x.Length;
double output = 0;
for (int i = 0; i < len; i++)
output += x[i];
return output / len;
}
public static double[] minus(double[] x, double[] y)
{ // subtract 2 vectors (element by element)
int len = x.Length;
if (len != y.Length)
{
if (y.Length == 1)
return minus(x,y[0]);
else
throw new ArgumentException("minus - input vectors must have the same length");
}
double[] output = new double[len];
for (int i = 0; i < len; i++)
{
output[i] = x[i] - y[i];
}
return output;
}
public static double[] minus(double[] x, double y)
{ // subtract a scalar from a vector
int len = x.Length;
double[] output = new double[len];
for (int i = 0; i < len; i++)
{
output[i] = x[i] - y;
}
return output;
}

Last edited by kevind; January 9th, 2012 at 11:27 PM .

The following 6 users say Thank You to kevind for this post:

April 10th, 2015, 09:22 AM
#12 (permalink )
Elite Member

Krakow Poland

Futures Experience: Intermediate

Platform: NinjaTrader

Favorite Futures: Stocks

Posts: 19 since Oct 2014

Thanks: 11 given,
1
received

kevind
Here's the ALGLIB file. It's a

NinjaTrader "Indicator", but it's just the entire ALGLIB 3.4 C# math library merged into a single file, put into a namespace called "ALGLIB". It's callable from strategies also.

I don't have time right now to give much of an example, but below are some static examples of matrix operations (least-squares fit, matrix multiplication, inverse, pseudoinverse, left divide & right divide, an OLS routine originally from MATLAB, and some helper routines). This (static) approach is easy to use (no instantiation required); it's fine for indicators but it falls down if you use it in a strategy which you run through the optimizer. That's because "new" is computationally expensive in C#, if you do "new" repeatedly it hits your performance and fills up all your memory. The alternative is to instantiate classes which allocate storage once and then re-use it, which helps a lot but it's more complicated and you still can eventually hit a wall for other reasons.

The code below is a bunch of copy & pastes; I think everything needed to compile it is here, but I haven't tested to be sure.

To use it:

Code
using ALGLIB;
public static void lsfitlinear(double[] y, double[,] x, out double[] beta)
{
int nobs = y.Length;
if (nobs != x.GetLength(0))
throw new ArgumentException("ols: y and x must have same length");
int nvar = x.GetLength(1);
int info;
alglib.lsfitreport rep;
alglib.lsfitlinear(y, x, nobs, nvar, out info, out beta, out rep);
}
public static double[,] mtimes(double[,] a, double[,] b, bool transposeA, bool transposeB)
{
int ra = a.GetLength(transposeA ? 1 : 0);
int ca = a.GetLength(transposeA ? 0 : 1);
int cb = b.GetLength(transposeB ? 0 : 1);
if (ca != b.GetLength(transposeB ? 1 : 0))
throw new ArgumentException("mtimes: columns(matrix a) must equal rows(matrix b)");
double [,] output = new double[ra,cb];
alglib.ablas.rmatrixgemm(ra,cb,ca,1,a,0,0,transposeA ? 1 : 0,b,0,0,transposeB ? 1 : 0,0,ref output,0,0);
return output;
}
public static double[,] inverse(double[,] a)
{
if (a.GetLength(0) != a.GetLength(1))
throw new ArgumentException("inverse: matrix must be square");
double[,] output = copy(a);
alglib.matinv.matinvreport rep = new alglib.matinv.matinvreport();
int info = new int();
alglib.matinv.rmatrixinverse(ref output,output.GetLength(0),ref info,rep);
return output;
}
public static double[,] copy(double[,] a)
{ // create a copy of a matrix
int r = a.GetLength(0);
int c = a.GetLength(1);
double[,] output = new double[r,c];
for (int i=0; i < r; i++)
{
for (int j=0; j < c; j++)
{
output[i,j] = a[i,j];
}
}
return output;
}
public static double[,] pseudoinverse(double[,] a)
{
int r = a.GetLength(0);
int c = a.GetLength(1);
int minrc = Math.Min(r,c);
double[] w = new double[minrc];
double[,] winv = new double[minrc,minrc];
double[,] u = new double[r,minrc];
double[,] vt = new double[minrc,c];
double[,] temp = new double[c,minrc];
double[,] output = new double[c,r];
alglib.svd.rmatrixsvd(a,r,c,1,1,2,ref w,ref u, ref vt);
// generate inverse winv of matrix W from the vector w from SVD
// to do this we transpose the matrix, and take the inverse of each vector element
// see: http://en.wikipedia.org/wiki/Singular_value_decomposition
// http://www.alglib.net/matrixops/general/svd.php
for (int i = 0; i < minrc; i++)
{
for (int j = 0; j < minrc; j++)
{
winv[i,j] = ((i == j) && (w[i] > alglib.math.machineepsilon)) ? 1.0/w[i] : 0.0; // inverse of w
}
}
alglib.ablas.rmatrixgemm(c,minrc,minrc,1,vt,0,0,1,winv,0,0,0,0,ref temp,0,0);
alglib.ablas.rmatrixgemm(c,r,minrc,1,temp,0,0,0,u,0,0,1,0,ref output,0,0);
return output;
}
public static double[,] mldivide(double[,] a, double[,] b)
{ // matrix left divide: inv(a)*b
int r = a.GetLength(0);
if (r != b.GetLength(0))
throw new ArgumentException ("mldivide: both matrices must have the same number of rows");
int ca = a.GetLength(1);
int cb = b.GetLength(1);
double[,] output = new double[ca,cb];
alglib.ablas.rmatrixgemm(ca,cb,r,1,pseudoinverse(a),0,0,0,b,0,0,0,0,ref output,0,0);
return output;
}
public static double[,] mrdivide(double[,] b, double[,] a)
{ // matrix right divide: b*inv(a)
int c = a.GetLength(1);
if (c != b.GetLength(1))
throw new ArgumentException ("mrdivide: both matrices must have the same number of columns");
int ra = a.GetLength(0);
int rb = b.GetLength(0);
double[,] output = new double[rb,ra];
alglib.ablas.rmatrixgemm(rb,ra,c,1,b,0,0,0,pseudoinverse(a),0,0,0,0,ref output,0,0);
return output;
}
public static void ols(double[] y, double[,] x, out double[] beta, out double[] yhat, out double[] residuals, out double sige, out double rsqr, out double rbar /*, out double[] tstats, out double[] bstd, out double dw*/)
{
// function results=ols(y,x)
// % PURPOSE: least-squares regression
// %---------------------------------------------------
// % USAGE: results = ols(y,x)
// % where: y = dependent variable vector (nobs x 1)
// % x = independent variables matrix (nobs x nvar)
// %---------------------------------------------------
// % RETURNS: a structure
// % results.meth = 'ols'
// % results.beta = bhat (nvar x 1)
// % results.tstat = t-stats (nvar x 1)
// % results.bstd = std deviations for bhat (nvar x 1)
// % results.yhat = yhat (nobs x 1)
// % results.resid = residuals (nobs x 1)
// % results.sige = e'*e/(n-k) scalar
// % results.rsqr = rsquared scalar
// % results.rbar = rbar-squared scalar
// % results.dw = Durbin-Watson Statistic
// % results.nobs = nobs
// % results.nvar = nvars
// % results.y = y data vector (nobs x 1)
// % results.bint = (nvar x2 ) vector with 95% confidence intervals on beta
// %---------------------------------------------------
int nobs = y.Length;
if (nobs != x.GetLength(0))
throw new ArgumentException("ols: y and x must have same length");
int nvar = x.GetLength(1);
int info;
yhat = new double[nobs];
alglib.lsfitreport rep;
alglib.lsfitlinear(y,x,nobs,nvar,out info, out beta, out rep); // least squares fit
alglib.ablas.rmatrixmv(x.GetLength(0),x.GetLength(1),x,0,0,0,beta,0,ref yhat,0); // calculate yhat (predicted values of Y)
residuals = minus(y,yhat);
double sigu = mtimes(residuals,residuals);
sige = sigu/System.Convert.ToDouble(nobs-nvar);
double[] ym = minus(y,mean(y));
double rsqr1 = sigu;
double rsqr2 = mtimes(ym,ym);
rsqr = 1.0 - rsqr1/rsqr2; // r-squared
rsqr1 = rsqr1/System.Convert.ToDouble(nobs-nvar);
rsqr2 = rsqr2/System.Convert.ToDouble(nobs-1);
rbar = rsqr2 == 0 ? rsqr : 1.0 - rsqr1/rsqr2;
// bstd, dw, tstat, bint not implemented
}
public static double mtimes(double[] x, double[] y)
{ // multiply 2 vectors (matrix multiply on 1-row and 1-col matrices)
int lx = x.Length;
if (lx != y.Length)
throw new ArgumentException("mtimes - input vectors must have the same length");
double output = 0;
for (int i = 0; i < lx; i++)
output += x[i] * y[i];
return output;
}
public static double mean(double[] x)
{ // mean of a vector
int len = x.Length;
double output = 0;
for (int i = 0; i < len; i++)
output += x[i];
return output / len;
}
public static double[] minus(double[] x, double[] y)
{ // subtract 2 vectors (element by element)
int len = x.Length;
if (len != y.Length)
{
if (y.Length == 1)
return minus(x,y[0]);
else
throw new ArgumentException("minus - input vectors must have the same length");
}
double[] output = new double[len];
for (int i = 0; i < len; i++)
{
output[i] = x[i] - y[i];
}
return output;
}
public static double[] minus(double[] x, double y)
{ // subtract a scalar from a vector
int len = x.Length;
double[] output = new double[len];
for (int i = 0; i < len; i++)
{
output[i] = x[i] - y;
}
return output;
}

Great stuff.
Do you have an example how I could create a forward moving OLS indicator, based on say 20 bars i.e calculates the coefficients on the last 20 bars, moves forward a bar and calculates for the next 20 and so on? I'll build on that. Thanks.

April 13th, 2015, 03:26 PM
#13 (permalink )
Elite Member

San Diego, CA

Futures Experience: Intermediate

Platform: x_trader pro, Excel

Broker/Data: Advantage

Favorite Futures: CL

Posts: 47 since Sep 2009

Thanks: 24 given,
24
received

Farmer George
Great stuff.
Do you have an example how I could create a forward moving OLS indicator, based on say 20 bars i.e calculates the coefficients on the last 20 bars, moves forward a bar and calculates for the next 20 and so on? I'll build on that. Thanks.

George, what platform / language are you using? Data-series are a native object to Ninjatrader .
If you are building your own app, the quickest way to go about this would be to include a time-series library in to your project.
If using C# or F#, check out Deedle
If using python, check out pandas
They both have examples of deriving rolling bar calculations.

April 13th, 2015, 03:54 PM
#14 (permalink )
Elite Member

Krakow Poland

Futures Experience: Intermediate

Platform: NinjaTrader

Favorite Futures: Stocks

Posts: 19 since Oct 2014

Thanks: 11 given,
1
received

baywolf
George, what platform / language are you using? Data-series are a native object to

Ninjatrader .

If you are building your own app, the quickest way to go about this would be to include a time-series library in to your project.

If using C# or F#, check out

Deedle
If using python, check out

pandas
They both have examples of deriving rolling bar calculations.

Baywolf, many thanks for that, appreciated.
I'm on fast learning curve with C#. I will use NinjaScript to create an indicator with a rolling Beta , regressing a share's returns on the 'market' returns. Also, I'll estimate and look at the hedge ratio between pairs of stocks/segments. I still think that there may be opportunities in pairs trading .
George

April 13th, 2015, 11:56 PM
#15 (permalink )
Trading Apprentice

Toronto, Canada

Futures Experience: Beginner

Platform: MarketDelta

Favorite Futures: Bitcoin

Posts: 6 since Apr 2015

Thanks: 0 given,
1
received

How about looking at the GNU/GSL library!!

April 14th, 2015, 07:11 AM
#16 (permalink )
Elite Member

Krakow Poland

Futures Experience: Intermediate

Platform: NinjaTrader

Favorite Futures: Stocks

Posts: 19 since Oct 2014

Thanks: 11 given,
1
received

markmycoin
How about looking at the GNU/GSL library!!

Ok, thanks. Had a quick look. Presumably one can get C# code?

Thread Tools
Search this Thread

Elite only
Elite only
Jan 18
Jan 23
Elite only
Elite only
Elite only

All times are GMT -4. The time now is 02:33 AM .