1 Introduction 1.1 What is ALGLIB 1.2 ALGLIB license 1.3 Documentation license 1.4 Reference Manual and User Guide 1.5 Acknowledgements 2 ALGLIB structure 2.1 Packages 2.2 Subpackages 2.3 Open Source and Commercial versions 3 Compatibility 4 Compiling ALGLIB 4.1 Adding to your project 5 Working with commercial version 5.1 Benefits of commercial version 5.2 Using multithreading 5.2.1 SMT (CMT/hyper-threading) issues 5.3 Using native computational core 5.3.1 Introduction 5.3.2 Linking with native core 5.3.3 HPC core vs. Managed core 6 Using ALGLIB 6.1 Thread-safety 6.2 ALGLIB classes 6.3 Datatypes 6.4 Functions 6.5 Using functions: 'expert' and 'friendly' interfaces 7 Advanced topics 7.1 Testing ALGLIB 7.2 Working with VERY large arrays 8 ALGLIB reference manual
ALGLIB is a cross-platform numerical analysis and data mining library. It supports several programming languages (C++, C#, Pascal, VBA) and several operating systems (Windows, *nix family).
ALGLIB features include:
ALGLIB Project (the company behind ALGLIB) delivers to you several editions of ALGLIB:
Free Edition is a serial version without multithreading support or extensive low-level optimizations (generic C or C# code). Commercial Edition is a heavily optimized version of ALGLIB. It supports multithreading, it was extensively optimized, and (on Intel platforms) - our commercial users may enjoy precompiled version of ALGLIB which internally calls Intel MKL to accelerate low-level tasks. We obtained license from Intel corp., which allows us to integrate Intel MKL into ALGLIB, so you don't have to buy separate license from Intel.
ALGLIB Free Edition is distributed under GPL 2+, GPL license version 2 or at your option any later version. A copy of the GNU General Public License is available at http://www.fsf.org/licensing/licenses
ALGLIB Commercial Edition is distributed under license which is friendly to commericial users. A copy of the commercial license can be found at http://www.alglib.net/commercial.php.
This reference manual is licensed under BSD-like documentation license:
Copyright 1994-2009 Sergey Bochkanov, ALGLIB Project. All rights reserved.
Redistribution and use of this document (ALGLIB Reference Manual) with or without modification, are permitted provided that such redistributions will retain the above copyright notice, this condition and the following disclaimer as the first (or last) lines of this file.

ALGLIB Project provides two sources of information: ALGLIB Reference Manual (this document) and ALGLIB User Guide.
ALGLIB Reference Manual contains full description of all publicly accessible ALGLIB units accompanied with examples.
Reference Manual is focused on the source code: it documents units, functions, structures and so on.
If you want to know what unit YYY
can do or what subroutines unit ZZZ
contains Reference Manual is a place to go.
Free software needs free documentation - that's why ALGLIB Reference Manual is licensed under BSD-like documentation license.
Additionally to the Reference Manual we provide you User Guide. User Guide is focused on more general questions: how fast ALGLIB is? how reliable it is? what are the strong and weak sides of the algorithms used? We aim to make ALGLIB User Guide an important source of information both about ALGLIB and numerical analysis algorithms in general. We want it to be a book about algorithms, not just software documentation. And we want it to be unique - that's why ALGLIB User Guide is distributed under less-permissive personal-use-only license.
ALGLIB was not possible without contribution of the next open source projects:
We also want to thank developers of the Intel's local development center (Nizhny Novgorod branch) for their help during MKL integration.
ALGLIB for C# is a pure C# 2.0 library automatically generated by code generation tools developed within ALGLIB project. Pre-3.0 versions of ALGLIB included more than 100 units, but it was difficult to work with such large number of files. Since ALGLIB 3.0 all units are merged into 11 packages and two support units:
One package may rely on other ones, but we have tried to reduce number of dependencies.
Every package relies on ap.cs
and many packages rely on alglibinternal.cs
.
But many packages require only these two to work, and many other packages need significantly less than 13 packages.
For example, statistics.cs
requires two packages mentioned above and only one additional package - specialfunctions.cs
.
There is one more concept to learn - subpackages.
Every package was created from several source files.
For example (as of ALGLIB 3.0.0), linalg.cs
was created by merging together 14 .cs files.
These files provide different functionality: one of them calculates triangular factorizations, another generates random matrices, and so on.
We've merged source code, but what to do with their documentation?
Of course, we can merge their documentation (as we've merged units) in one big list of functions and data structures, but such list will be hard to read. Instead, we have decided to merge source code, but leave documentation separate.
If you look at the list of ALGLIB packages, you will see that each package includes several subpackages.
For example, linalg.cs
includes trfac
, svd
, evd
and other subpackages.
These subpackages do no exist as separate files, namespaces or other entities.
They are just subsets of one large file, one large class, which provide significantly different functionality.
They have separate documentation sections, but if you want to use svd
subpackage, you have to compile linalg.cs
, not svd.cs
.
ALGLIB comes in two versions - open source (GPL-licensed) and commercial (closed source) one. Both versions have same functionality, i.e. may solve same set of problems. However, commercial version differs from open source one in following aspects:
This documentation applies to both versions of ALGLIB. Parts specific to commercial versions are explicitly marked as such.
ALGLIB for C# is compatible with:
Adding ALGLIB to your project is easy - just pick packages you need and... add them to your project! It will work without any additional settings.
As you see, ALGLIB has no project files. Why? There are several reasons:
In any case, compiling ALGLIB is so simple that even without project file you can do it in several minutes.
Commercial version of ALGLIB for C# features three important improvements over open source one:
Commercial version of ALGLIB includes out-of-the-box support for multithreading. Many (not all) computationally intensive problems can be solved in multithreaded mode. You should read comments on specific ALGLIB functions to determine what can be multithreaded and what can not.
ALGLIB does not depend on vendor/compiler support for technologies like OpenMP/MPI/... ALGLIB uses managed threads provided by NET framework. It gives ALGLIB unprecedented portability across operating systems and compilers. ALGLIB does not depend on presence of any custom multithreading library or compiler support for any multithreading technology.
If you want to use multithreaded capabilities of ALGLIB, you should:
Let explain it in more details...
1.
ALGLIB automatically determines number of cores on application startup.
By default, ALGLIB uses all available cores except for one.
Say, on 4-core system it will use three cores - unless being told to use more or less.
It will keep your system responsive during lengthy computations.
Such behavior may be changed with setnworkers()
call:
alglib.setnworkers(0)
= use all coresalglib.setnworkers(-1)
= leave one core unusedalglib.setnworkers(-2)
= leave two cores unusedalglib.setnworkers(+2)
= use 2 cores (even if you have more)
You may want to specify maximum number of worker threads during compile time
by means of preprocessor definition AE_NWORKERSX
(say, AE_NWORKERS2, AE_NWORKERS4, AE_NWORKERS5).
You can add this definition to compiler command line or change corresponding project settings in your IDE.
Here X can be any positive number.
ALGLIB will use exactly X worker threads, unless being told to use less by setnworkers()
call.
2. When you use commercial edition of ALGLIB, you may choose between serial and multithreaded versions of SMP-capable functions:
smp_
prefix in its name) creates (or wakes up) worker threads,
inserts task in the worker queue, and waits for completion of the task.
All processing is done in context of worker thread(s).You should carefully decide what version of function to use. Starting/stopping worker thread costs tens of thousands of CPU cycles. Thus you won't get multithreading speedup on small computational problems.
Simultaneous multithreading (SMT) also known as Hyper-threading (Intel) and Cluster-based Multithreading (AMD) is a CPU design where several (usually two) logical cores share resources of one physical core. Say, on dual-core system with 2x HT scale factor you will see 4 logical cores. Each pair of these 4 cores, however, share same hardware resources. Thus, you may get only marginal speedup when running highly optimized software which fully utilizes CPU resources.
Say, if one thread occupies floating-point unit, another thread on the same physical core may work with integer numbers at the same time without any performance penalties. In this case you may get some speedup due to having additional cores. But if both threads keep FPU unit 100% busy, they won't get any multithreaded speedup.
So, if 2 math-intensive threads are dispatched by OS scheduler to different physical cores, you will get 2x speedup due to use of multithreading. But if these threads are dispatched to different logical cores - but same physical core - you won't get any speedup at all! One physical core will be 100% busy, and another one will be 100% idle. From the other side, if you start four threads instead of two, your system will be 100% utilized independently of thread scheduling details.
Let we stress it one more time - multithreading speedup on SMT systems is highly dependent on number of threads you are running and decisions made by OS scheduler. It is not 100% deterministic! With "true SMP" when you run 2 threads, you get 2x speedup (or 1.95, or 1.80 - it depends on algorithm, but this factor is always same). With SMT when you run 2 threads you may get your 2x speedup - or no speedup at all. Modern OS schedulers do a good job on single-socket hardware, but even in this "simple" case they give no guarantees of fair distribution of hardware resources. And things become a bit tricky when you work with multi-socket hardware. On SMT systems the only guaranteed way to 100% utilize your CPU is to create as many worker threads as there are logical cores. In this case OS scheduler has no chance to make its work in a wrong way.
Commercial version of ALGLIB for C# includes high-performance native computational core. It supports multithreading and it was linked with Intel MKL - highly optimized linear algebra and FFT library from Intel. Many ALGLIB functions are accelerated by MKL.
Native core provides same functionality as managed one and has 100% identical interface. It allows you to switch between native and managed cores as you wish. You may choose portability of fully managed code - or high performance of native core. You may change cores during application development - without ever changing code which calls ALGLIB!
As of ALGLIB 3.8.1, HPC core can work only on Windows platforms (both x86 and x86_64). However, we plan to introduce *nix versions of HPC core in future ALGLIB releases!
Native computational core is located in hpc-windows
subdirectory of the ALGLIB distribution.
This directory contains following files:
alglib32_hpc.dll
- 32-bit native DLL which includes
C version of ALGLIB compiled with best optimization settings and linked with Intel MKL
(MKL is statically linked and included in resulting DLL)
alglib64_hpc.dll
- 64-bit native DLL with same functionality
alglib_hpc.cs
- C# wrapper file which automatically detects platform we are running on (x86 or x64)
and loads appropriate DLL (32-bit or 64-bit one).
First step is to add C# wrapper file alglib_hpc.cs
to your project.
Just compile it with the rest of your application as one more unit.
This wrapper provides exactly same set of classes and methods as managed version of ALGLIB.
However, all method calls are silently redirected to native DLL,
and all ALGLIB objects now wrappers around native data structures.
IMPORTANT: you will have to allow "unsafe" code in your application in order to compile ALGLIB wrapper layer.
In can be done by specifying "/insafe+"
at command line -
or by changing appropriate field in your project settings.
ALGLIB native DLL is automatically loaded during first call of any ALGLIB function. Wrapper layer detects platform it is running on (x86 or x64) and chooses appropriate DLL (32-bit or 64-bit one). Thus your application can be easily moved between 32-bit and 64-bit hosts without recompilation.
ALGLIB native DLL's are automatically loaded during first call of ALGLIB functions,
but you may enforce this process by calling alglib.activatealglibcore
method.
It can be useful if you running some kind of benchmark -
or just want to perform all "heavy" initialization during application startup.
Because HPC ALGLIB is a wrapper around native DLL's, these DLL's must be always accompany your application.
Wrapper layer tries to find native DLL's in the application's home directory.
It works well for "standalone" executables,
but this method may fail when wrapper is a part of assembly which was placed to the GAC.
In such cases you may manually specify path to ALGLIB native DLL's by calling alglib.addsearchpath
method.
This method must be called prior to first call of ALGLIB functions.
After all steps above are done you may start using commercial version of ALGLIB!
It was told before that HPC core provides same functionality as its managed counterpart. It means that any application developed with managed core will work with HPC core. However, there are several notable difference between two cores. In this section we will outline these differences.
Performance. HPC core has better performance on large-scale problems (2x-4x improvement for most tasks), but it is less suited for small-scale computational problems. Say, calling HPC linear algebra functions involves copying of managed input arrays to unmanaged memory - and similar extraction of result from unmanaged memory back to managed world. Performance impact is negligible when you work with 1000x1000 matrices (multiplying them takes much longer then copying them), but if you want to multiply 3x3 matrices - you will get better results with managed core, which performs all work in managed memory.
Reproducibility of numerical results. When your NET application works with floating point data, it may produce a bit different results on different computers. It is normal situation which is caused by presence of rounding errors, which may have different magnitudes on different machines. Same situation can be observed with HPC and managed ALGLIB cores - same algorithm may produce a bit different results in different computational cores. Again, it is not algorithmic error, it just consequence of the fact that CPU works with finite-precision numbers.
Memory consumption and management. Consider situation when you create ALGLIB object which internally allocates large chunk of ALGLIB-controlled memory, and then you destroy this object.
Working with large data
alglib.kdtree tree; alglib.kdtreebuild(onegigabyteofdata, ..., out tree); // Pass 1GB of data to object tree = null; // Destroy object
What will do managed and native versions of ALGLIB after last line? Pure NET version of ALGLIB stores all its data in managed memory and it will handle this situation just as normal NET application. Gigabyte-sized memory chunk will be forgotten and reclaimed later during garbage collection. What is most important, NET framework will be aware that we have 1GB of memory which waits to be collected.
However, situation is different with HPC ALGLIB.
This version of ALGLIB stores data in the unmanaged memory.
In HPC ALGLIB, kdtree
object stores just opaque pointer to unmanaged data structures.
One gigabyte of data is copied to unmanaged memory, where it is controlled by ALGLIB -
and NET garbage collector is no longer aware about this 1GB-sized block.
It knows about ALGLIB object which is responsible for this data block,
but from the GC's point of view, size of this object is 8 bytes (size of pointer).
It has no way of knowing about large data structures hidden behind this pointer.
Of course, this memory won't leak.
After last line kdtree
object will become candidate for the next garbage collection run.
It has destructor which will be called during garbage collection,
so ALGLIB will be able to free all memory it is responsible for.
The problem is that you don't know when it will be done.
Because GC do not know about actual size of ALGLIB object,
it may decide that we do not have to perform garbage collection right now.
So it may be beneficial to enforce GC run, as demonstrated below:
Working with large data : solution #1
alglib.kdtree tree; alglib.kdtreebuild(onegigabyteofdata, ..., out tree); // Pass 1GB of data to object tree = null; // Destroy object System.GC.Collect(); // Enforce immediate GC
Such solution works in most cases.
However, because GC is unaware about actual size of tree
object,
sometimes it may decide to put such tiny object (8 bytes!) on hold until the next run.
And next.
And next.
Official documentation from Microsoft says that GC.Collect
collects all objects,
but in practice some objects can survive GC run.
To handle such situations, HPC ALGLIB includes special method: alglib.deallocateimmediately
,
which enforces immediate deallocation of internal data structures owned by ALGLIB object.
After method was called on ALGLIB object, all its memory is released and object becomes unusable.
In managed edition of ALGLIB this method does nothing, just sets variable it accepts to null.
Working with large data : solution #2
alglib.kdtree tree; alglib.kdtreebuild(onegigabyteofdata, ..., out tree); // Pass 1GB of data to object alglib.deallocateimmediately(ref tree); // deallocated and set to null
Both open source and commercial versions of ALGLIB are 100% thread-safe as long as different user threads work with different instances of objects/arrays. Thread-safety is guaranteed by having no global shared variables.
However, any kind of sharing ALGLIB objects/arrays between different threads is potentially hazardous. Even when this object is seemingly used in read-only mode!
Say, you use ALGLIB neural network NET to process two input vectors X0 and X1,
and get two output vectors Y0 and Y1.
You may decide that neural network is used in read-only mode which does not change state of NET,
because output is written to distinct arrays Y.
Thus, you may want to process these vectors from parallel threads.
But it is not read-only operation, even if it looks like this!
Neural network object NET allocates internal temporary buffers, which are modified by neural processing functions.
Thus, sharing one instance of neural network between two threads is thread-unsafe!
All ALGLIB functionality is provided by one large class - alglib
.
ALGLIB defines several "basic" datatypes (types which are used across many packages) and many package-specific datatypes.
"Basic" datatypes are defined in ap.cs
.
Here is short list:
alglib.complex
- double precision complex datatype, which supports arithmetic operationsalglib.alglibexception
- exception which is thrown by libraryalglib.ndimensional_func
- delegate used by optimizers/solversalglib.ndimensional_grad
- delegate used by optimizers/solversalglib.ndimensional_hess
- delegate used by optimizers/solversalglib.ndimensional_jac
- delegate used by optimizers/solversalglib.ndimensional_pfunc
- delegate used by optimizers/solversalglib.ndimensional_pgrad
- delegate used by optimizers/solversalglib.ndimensional_phess
- delegate used by optimizers/solversalglib.ndimensional_rep
- delegate used by optimizers/solversalglib.ndimensional_ode_rp
- delegate used by ODE solveralglib.integrator1_func
- delegate used by 1-dimensional integratorPackage-specific datatypes are classes which can be divided into two distinct groups:
The most important "basic" functions from alglib
class are:
alglib.randomreal()
- returns random real number from [0,1)alglib.randominteger(mx)
- returns random integer number from [0,nx)alglib.sqr(x)
- returns square of real xalglib.csqr(z)
- returns square of complex zalglib.conj(z)
- returns conjugate of complex zalglib.abscomplex(z)
- returns magnitude of complex zMost ALGLIB functions provide two interfaces: 'expert' and 'friendly'. What is the difference between two? When you use 'friendly' interface, ALGLIB:
When you use 'expert' interface, ALGLIB:
Here are several examples of 'friendly' and 'expert' interfaces:
double[] x = new double[]{0,1,2,3}; double[] y = new double[]{1,5,3,9}; double[] y2 = new double[]{1,5,3,9,0}; alglib.spline1dinterpolant s; alglib.spline1dbuildlinear(x, y, 4, out s); // 'expert' interface is used alglib.spline1dbuildlinear(x, y, out s); // 'friendly' interface - input size is // automatically determined alglib.spline1dbuildlinear(x, y2, 4, out s); // y2.length() is 5, but it will work alglib.spline1dbuildlinear(x, y2, out s); // it won't work because sizes of x and y2 // are inconsistent
'Friendly' interface - matrix semantics:
double[,] a; alglib.matinvreport rep; alglib.ae_int_t info; // // 'Friendly' interface: spdmatrixinverse() accepts and returns symmetric matrix // // symmetric positive definite matrix a = new double[,]{{2,1},{1,2}}; // after this line A will contain [[0.66,-0.33],[-0.33,0.66]] // which is symmetric too alglib.spdmatrixinverse(ref a, out info, out rep); // you may try to pass nonsymmetric matrix a = new double[,]{{2,1},{0,2}}; // but exception will be thrown in such case alglib.spdmatrixinverse(ref a, out info, out rep);
Same function but with 'expert' interface:
double[,] a; alglib.matinvreport rep; alglib.ae_int_t info; // // 'Expert' interface, spdmatrixinverse() // // only upper triangle is used; a[1][0] is initialized by NAN, // but it can be arbitrary number a = new double[,]{{2,1},{Double.NaN,2}}; // after this line A will contain [[0.66,-0.33],[NAN,0.66]] // only upper triangle is modified alglib.spdmatrixinverse(ref a, 2 /* N */, true /* upper triangle is used */, out info, out rep);
There are two test suites in ALGLIB: computational tests and interface tests.
Computational tests are located in /tests/test_c.cs
.
They are focused on numerical properties of algorithms, stress testing and "deep" tests (large automatically generated problems).
They require significant amount of time to finish (tens of minutes).
Interface tests are located in /tests/test_i.cs
.
These tests are focused on ability to correctly pass data between computational core and caller, ability to detect simple problems in inputs,
and on ability to at least compile ALGLIB with your compiler.
They are very fast (about a minute to finish including compilation time).
Running test suite is easy - just
test_c.cs
or test_i.cs
) along with the rest of the libraryHowever, there is no strong reasons to launch tests before using ALGLIB. It was tested with many different compiler settings, and the very nature of .NET framework and its inherent portability allows us to say that ALGLIB for .NET will work on any system.
Sometimes ALGLIB users have to work with very large computational problems. For example, RBF interpolation problem with 400.000 points may require several gigabytes of memory for temporary variables (mostly for sparse matrix of coefficients). The problem is that .NET framework has limitations on array size even on 64-bit systems: 1) array size must be less than 2GB, 2) maximum index in any single dimension is about 2^31, 3) maximum number of elements in array is UInt32MaxValue.
The most important restriction is (1) - 2GB limit on array size. Say, for double precision values we can not have array with more than 2^28 (about 268.000.000) elements, which is 8 times smaller than limitation imposed by (2) and (3).
However, NET 4.5 has ability to allocate arrays larger than 2GB (with all other restrictions being present).
In order to do so you have to use gcAllowVeryLargeObjects=true
in your configuration schema.
With this setting turned on you can allocate double precision arrays with up to 2^31 elements,
which use up to 16GB of memory.
AlglibMisc package | ||
| hqrnd | High quality random numbers generator | |
| nearestneighbor | Nearest neighbor search: approximate and exact | |
| xdebug | Debug functions to test ALGLIB interface generator | |
DataAnalysis package | ||
| bdss | Basic dataset functions | |
| clustering | Clustering functions (hierarchical, k-means++) | |
| datacomp | Backward compatibility functions | |
| dforest | Decision forest classifier (regression model) | |
| filters | Different filters used in data analysis | |
| lda | Linear discriminant analysis | |
| linreg | Linear models | |
| logit | Logit models | |
| mcpd | Markov Chains for Population/proportional Data | |
| mlpbase | Basic functions for neural networks | |
| mlpe | Basic functions for neural ensemble models | |
| mlptrain | Neural network training | |
| pca | Principal component analysis | |
DiffEquations package | ||
| odesolver | Ordinary differential equation solver | |
FastTransforms package | ||
| conv | Fast real/complex convolution | |
| corr | Fast real/complex cross-correlation | |
| fft | Real/complex FFT | |
| fht | Real Fast Hartley Transform | |
Integration package | ||
| autogk | Adaptive 1-dimensional integration | |
| gkq | Gauss-Kronrod quadrature generator | |
| gq | Gaussian quadrature generator | |
Interpolation package | ||
| idwint | Inverse distance weighting: interpolation/fitting | |
| lsfit | Linear and nonlinear least-squares solvers | |
| polint | Polynomial interpolation/fitting | |
| pspline | Parametric spline interpolation | |
| ratint | Rational interpolation/fitting | |
| rbf | Scattered 2/3-dimensional interpolation with RBF models | |
| spline1d | 1D spline interpolation/fitting | |
| spline2d | 2D spline interpolation | |
| spline3d | 3D spline interpolation | |
LinAlg package | ||
| ablas | Level 2 and Level 3 BLAS operations | |
| bdsvd | Bidiagonal SVD | |
| evd | Eigensolvers | |
| inverseupdate | Sherman-Morrison update of the inverse matrix | |
| matdet | Determinant calculation | |
| matgen | Random matrix generation | |
| matinv | Matrix inverse | |
| normestimator | Estimates norm of the sparse matrix (from below) | |
| ortfac | Real/complex QR/LQ, bi(tri)diagonal, Hessenberg decompositions | |
| rcond | Condition number estimate | |
| schur | Schur decomposition | |
| sparse | Sparse matrices | |
| spdgevd | Generalized symmetric eigensolver | |
| svd | Singular value decomposition | |
| trfac | LU and Cholesky decompositions | |
Optimization package | ||
| minbleic | Bound constrained optimizer with additional linear equality/inequality constraints | |
| mincg | Conjugate gradient optimizer | |
| mincomp | Backward compatibility functions | |
| minlbfgs | Limited memory BFGS optimizer | |
| minlm | Improved Levenberg-Marquardt optimizer | |
| minqp | Quadratic programming with bound and linear equality/inequality constraints | |
Solvers package | ||
| densesolver | Dense linear system solver | |
| lincg | Sparse linear CG solver | |
| linlsqr | Sparse linear LSQR solver | |
| nleq | Solvers for nonlinear equations | |
SpecialFunctions package | ||
| airyf | Airy functions | |
| bessel | Bessel functions | |
| betaf | Beta function | |
| binomialdistr | Binomial distribution | |
| chebyshev | Chebyshev polynomials | |
| chisquaredistr | Chi-Square distribution | |
| dawson | Dawson integral | |
| elliptic | Elliptic integrals | |
| expintegrals | Exponential integrals | |
| fdistr | F-distribution | |
| fresnel | Fresnel integrals | |
| gammafunc | Gamma function | |
| hermite | Hermite polynomials | |
| ibetaf | Incomplete beta function | |
| igammaf | Incomplete gamma function | |
| jacobianelliptic | Jacobian elliptic functions | |
| laguerre | Laguerre polynomials | |
| legendre | Legendre polynomials | |
| normaldistr | Normal distribution | |
| poissondistr | Poisson distribution | |
| psif | Psi function | |
| studenttdistr | Student's t-distribution | |
| trigintegrals | Trigonometric integrals | |
Statistics package | ||
| basestat | Mean, variance, covariance, correlation, etc. | |
| correlationtests | Hypothesis testing: correlation tests | |
| jarquebera | Hypothesis testing: Jarque-Bera test | |
| mannwhitneyu | Hypothesis testing: Mann-Whitney-U test | |
| stest | Hypothesis testing: sign test | |
| studentttests | Hypothesis testing: Student's t-test | |
| variancetests | Hypothesis testing: F-test and one-sample variance test | |
| wsr | Hypothesis testing: Wilcoxon signed rank test | |
ablas
subpackage| ablas_d_gemm | Matrix multiplication (single-threaded) | |
| ablas_d_syrk | Symmetric rank-K update (single-threaded) | |
| ablas_smp_gemm | Matrix multiplication (multithreaded) | |
| ablas_smp_syrk | Symmetric rank-K update (multithreaded) |
cmatrixcopy
function/************************************************************************* Copy Input parameters: M - number of rows N - number of columns A - source matrix, MxN submatrix is copied and transposed IA - submatrix offset (row index) JA - submatrix offset (column index) B - destination matrix, must be large enough to store result IB - submatrix offset (row index) JB - submatrix offset (column index) *************************************************************************/public static void alglib.cmatrixcopy( int m, int n, complex[,] a, int ia, int ja, ref complex[,] b, int ib, int jb)
cmatrixgemm
function/************************************************************************* *************************************************************************/public static void alglib.cmatrixgemm( int m, int n, int k, complex alpha, complex[,] a, int ia, int ja, int optypea, complex[,] b, int ib, int jb, int optypeb, complex beta, ref complex[,] c, int ic, int jc) public static void alglib.smp_cmatrixgemm( int m, int n, int k, complex alpha, complex[,] a, int ia, int ja, int optypea, complex[,] b, int ib, int jb, int optypeb, complex beta, ref complex[,] c, int ic, int jc)
cmatrixlefttrsm
function/************************************************************************* *************************************************************************/public static void alglib.cmatrixlefttrsm( int m, int n, complex[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref complex[,] x, int i2, int j2) public static void alglib.smp_cmatrixlefttrsm( int m, int n, complex[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref complex[,] x, int i2, int j2)
cmatrixmv
function/************************************************************************* Matrix-vector product: y := op(A)*x INPUT PARAMETERS: M - number of rows of op(A) M>=0 N - number of columns of op(A) N>=0 A - target matrix IA - submatrix offset (row index) JA - submatrix offset (column index) OpA - operation type: * OpA=0 => op(A) = A * OpA=1 => op(A) = A^T * OpA=2 => op(A) = A^H X - input vector IX - subvector offset IY - subvector offset Y - preallocated matrix, must be large enough to store result OUTPUT PARAMETERS: Y - vector which stores result if M=0, then subroutine does nothing. if N=0, Y is filled by zeros. -- ALGLIB routine -- 28.01.2010 Bochkanov Sergey *************************************************************************/public static void alglib.cmatrixmv( int m, int n, complex[,] a, int ia, int ja, int opa, complex[] x, int ix, ref complex[] y, int iy)
cmatrixrank1
function/************************************************************************* Rank-1 correction: A := A + u*v' INPUT PARAMETERS: M - number of rows N - number of columns A - target matrix, MxN submatrix is updated IA - submatrix offset (row index) JA - submatrix offset (column index) U - vector #1 IU - subvector offset V - vector #2 IV - subvector offset *************************************************************************/public static void alglib.cmatrixrank1( int m, int n, ref complex[,] a, int ia, int ja, ref complex[] u, int iu, ref complex[] v, int iv)
cmatrixrighttrsm
function/************************************************************************* *************************************************************************/public static void alglib.cmatrixrighttrsm( int m, int n, complex[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref complex[,] x, int i2, int j2) public static void alglib.smp_cmatrixrighttrsm( int m, int n, complex[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref complex[,] x, int i2, int j2)
cmatrixsyrk
function/************************************************************************* *************************************************************************/public static void alglib.cmatrixsyrk( int n, int k, double alpha, complex[,] a, int ia, int ja, int optypea, double beta, ref complex[,] c, int ic, int jc, bool isupper) public static void alglib.smp_cmatrixsyrk( int n, int k, double alpha, complex[,] a, int ia, int ja, int optypea, double beta, ref complex[,] c, int ic, int jc, bool isupper)
cmatrixtranspose
function/************************************************************************* Cache-oblivous complex "copy-and-transpose" Input parameters: M - number of rows N - number of columns A - source matrix, MxN submatrix is copied and transposed IA - submatrix offset (row index) JA - submatrix offset (column index) B - destination matrix, must be large enough to store result IB - submatrix offset (row index) JB - submatrix offset (column index) *************************************************************************/public static void alglib.cmatrixtranspose( int m, int n, complex[,] a, int ia, int ja, ref complex[,] b, int ib, int jb)
rmatrixcopy
function/************************************************************************* Copy Input parameters: M - number of rows N - number of columns A - source matrix, MxN submatrix is copied and transposed IA - submatrix offset (row index) JA - submatrix offset (column index) B - destination matrix, must be large enough to store result IB - submatrix offset (row index) JB - submatrix offset (column index) *************************************************************************/public static void alglib.rmatrixcopy( int m, int n, double[,] a, int ia, int ja, ref double[,] b, int ib, int jb)
rmatrixenforcesymmetricity
function/************************************************************************* This code enforces symmetricy of the matrix by copying Upper part to lower one (or vice versa). INPUT PARAMETERS: A - matrix N - number of rows/columns IsUpper - whether we want to copy upper triangle to lower one (True) or vice versa (False). *************************************************************************/public static void alglib.rmatrixenforcesymmetricity( ref double[,] a, int n, bool isupper)
rmatrixgemm
function/************************************************************************* *************************************************************************/public static void alglib.rmatrixgemm( int m, int n, int k, double alpha, double[,] a, int ia, int ja, int optypea, double[,] b, int ib, int jb, int optypeb, double beta, ref double[,] c, int ic, int jc) public static void alglib.smp_rmatrixgemm( int m, int n, int k, double alpha, double[,] a, int ia, int ja, int optypea, double[,] b, int ib, int jb, int optypeb, double beta, ref double[,] c, int ic, int jc)
rmatrixlefttrsm
function/************************************************************************* *************************************************************************/public static void alglib.rmatrixlefttrsm( int m, int n, double[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref double[,] x, int i2, int j2) public static void alglib.smp_rmatrixlefttrsm( int m, int n, double[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref double[,] x, int i2, int j2)
rmatrixmv
function/************************************************************************* Matrix-vector product: y := op(A)*x INPUT PARAMETERS: M - number of rows of op(A) N - number of columns of op(A) A - target matrix IA - submatrix offset (row index) JA - submatrix offset (column index) OpA - operation type: * OpA=0 => op(A) = A * OpA=1 => op(A) = A^T X - input vector IX - subvector offset IY - subvector offset Y - preallocated matrix, must be large enough to store result OUTPUT PARAMETERS: Y - vector which stores result if M=0, then subroutine does nothing. if N=0, Y is filled by zeros. -- ALGLIB routine -- 28.01.2010 Bochkanov Sergey *************************************************************************/public static void alglib.rmatrixmv( int m, int n, double[,] a, int ia, int ja, int opa, double[] x, int ix, ref double[] y, int iy)
rmatrixrank1
function/************************************************************************* Rank-1 correction: A := A + u*v' INPUT PARAMETERS: M - number of rows N - number of columns A - target matrix, MxN submatrix is updated IA - submatrix offset (row index) JA - submatrix offset (column index) U - vector #1 IU - subvector offset V - vector #2 IV - subvector offset *************************************************************************/public static void alglib.rmatrixrank1( int m, int n, ref double[,] a, int ia, int ja, ref double[] u, int iu, ref double[] v, int iv)
rmatrixrighttrsm
function/************************************************************************* *************************************************************************/public static void alglib.rmatrixrighttrsm( int m, int n, double[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref double[,] x, int i2, int j2) public static void alglib.smp_rmatrixrighttrsm( int m, int n, double[,] a, int i1, int j1, bool isupper, bool isunit, int optype, ref double[,] x, int i2, int j2)
rmatrixsyrk
function/************************************************************************* *************************************************************************/public static void alglib.rmatrixsyrk( int n, int k, double alpha, double[,] a, int ia, int ja, int optypea, double beta, ref double[,] c, int ic, int jc, bool isupper) public static void alglib.smp_rmatrixsyrk( int n, int k, double alpha, double[,] a, int ia, int ja, int optypea, double beta, ref double[,] c, int ic, int jc, bool isupper)
rmatrixtranspose
function/************************************************************************* Cache-oblivous real "copy-and-transpose" Input parameters: M - number of rows N - number of columns A - source matrix, MxN submatrix is copied and transposed IA - submatrix offset (row index) JA - submatrix offset (column index) B - destination matrix, must be large enough to store result IB - submatrix offset (row index) JB - submatrix offset (column index) *************************************************************************/public static void alglib.rmatrixtranspose( int m, int n, double[,] a, int ia, int ja, ref double[,] b, int ib, int jb)
public static int Main(string[] args) { double[,] a = new double[,]{{2,1},{1,3}}; double[,] b = new double[,]{{2,1},{0,1}}; double[,] c = new double[,]{{0,0},{0,0}}; // // rmatrixgemm() function allows us to calculate matrix product C:=A*B or // to perform more general operation, C:=alpha*op1(A)*op2(B)+beta*C, // where A, B, C are rectangular matrices, op(X) can be X or X^T, // alpha and beta are scalars. // // This function: // * can apply transposition and/or multiplication by scalar to operands // * can use arbitrary part of matrices A/B (given by submatrix offset) // * can store result into arbitrary part of C // * for performance reasons requires C to be preallocated // // Parameters of this function are: // * M, N, K - sizes of op1(A) (which is MxK), op2(B) (which // is KxN) and C (which is MxN) // * Alpha - coefficient before A*B // * A, IA, JA - matrix A and offset of the submatrix // * OpTypeA - transformation type: // 0 - no transformation // 1 - transposition // * B, IB, JB - matrix B and offset of the submatrix // * OpTypeB - transformation type: // 0 - no transformation // 1 - transposition // * Beta - coefficient before C // * C, IC, JC - preallocated matrix C and offset of the submatrix // // Below we perform simple product C:=A*B (alpha=1, beta=0) // // IMPORTANT: this function works with preallocated C, which must be large // enough to store multiplication result. // int m = 2; int n = 2; int k = 2; double alpha = 1.0; int ia = 0; int ja = 0; int optypea = 0; int ib = 0; int jb = 0; int optypeb = 0; double beta = 0.0; int ic = 0; int jc = 0; alglib.rmatrixgemm(m, n, k, alpha, a, ia, ja, optypea, b, ib, jb, optypeb, beta, ref c, ic, jc); System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[4,3],[2,4]] // // Now we try to apply some simple transformation to operands: C:=A*B^T // optypeb = 1; alglib.rmatrixgemm(m, n, k, alpha, a, ia, ja, optypea, b, ib, jb, optypeb, beta, ref c, ic, jc); System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[5,1],[5,3]] System.Console.ReadLine(); return 0; }
public static int Main(string[] args) { // // rmatrixsyrk() function allows us to calculate symmetric rank-K update // C := beta*C + alpha*A'*A, where C is square N*N matrix, A is square K*N // matrix, alpha and beta are scalars. It is also possible to update by // adding A*A' instead of A'*A. // // Parameters of this function are: // * N, K - matrix size // * Alpha - coefficient before A // * A, IA, JA - matrix and submatrix offsets // * OpTypeA - multiplication type: // * 0 - A*A^T is calculated // * 2 - A^T*A is calculated // * Beta - coefficient before C // * C, IC, JC - preallocated input/output matrix and submatrix offsets // * IsUpper - whether upper or lower triangle of C is updated; // this function updates only one half of C, leaving // other half unchanged (not referenced at all). // // Below we will show how to calculate simple product C:=A'*A // // NOTE: beta=0 and we do not use previous value of C, but still it // MUST be preallocated. // int n = 2; int k = 1; double alpha = 1.0; int ia = 0; int ja = 0; int optypea = 2; double beta = 0.0; int ic = 0; int jc = 0; bool isupper = true; double[,] a = new double[,]{{1,2}}; // preallocate space to store result double[,] c = new double[,]{{0,0},{0,0}}; // calculate product, store result into upper part of c alglib.rmatrixsyrk(n, k, alpha, a, ia, ja, optypea, beta, ref c, ic, jc, isupper); // output result. // IMPORTANT: lower triangle of C was NOT updated! System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[1,2],[0,4]] System.Console.ReadLine(); return 0; }
public static int Main(string[] args) { // // In this example we assume that you already know how to work with // rmatrixgemm() function. Below we concentrate on its multithreading // capabilities. // // SMP edition of ALGLIB includes smp_rmatrixgemm() - multithreaded // version of rmatrixgemm() function. In the basic edition of ALGLIB // (GPL edition or commercial version without SMP support) this function // just calls single-threaded stub. So, you may call this function from // ANY edition of ALGLIB, but only in SMP edition it will work in really // multithreaded mode. // // In order to use multithreading, you have to: // 1) Install SMP edition of ALGLIB. // 2) This step is specific for C++ users: you should activate OS-specific // capabilities of ALGLIB by defining AE_OS=AE_POSIX (for *nix systems) // or AE_OS=AE_WINDOWS (for Windows systems). // C# users do not have to perform this step because C# programs are // portable across different systems without OS-specific tuning. // 3) Allow ALGLIB to know about number of worker threads to use: // a) autodetection (C++, C#): // ALGLIB will automatically determine number of CPU cores and // (by default) will use all cores except for one. Say, on 4-core // system it will use three cores - unless you manually told it // to use more or less. It will keep your system responsive during // lengthy computations. // Such behavior may be changed with setnworkers() call: // * alglib::setnworkers(0) = use all cores // * alglib::setnworkers(-1) = leave one core unused // * alglib::setnworkers(-2) = leave two cores unused // * alglib::setnworkers(+2) = use 2 cores (even if you have more) // b) manual specification (C++, C#): // You may want to specify maximum number of worker threads during // compile time by means of preprocessor definition AE_NWORKERS. // For C++ it will be "AE_NWORKERS=X" where X can be any positive number. // For C# it is "AE_NWORKERSX", where X should be replaced by number of // workers (AE_NWORKERS2, AE_NWORKERS3, AE_NWORKERS4, ...). // You can add this definition to compiler command line or change // corresponding project settings in your IDE. // // After you installed and configured SMP edition of ALGLIB, you may choose // between serial and multithreaded versions of SMP-capable functions: // * serial version works as usual, in the context of the calling thread // * multithreaded version (with "smp_" prefix) creates (or wakes up) worker // threads, inserts task in the worker queue, and waits for completion of // the task. All processing is done in context of worker thread(s). // // NOTE: because starting/stopping worker threads costs thousands of CPU cycles, // you should not use multithreading for lightweight computational problems. // // NOTE: some old POSIX-compatible operating systems do not support // sysconf(_SC_NPROCESSORS_ONLN) system call which is required in order // to automatically determine number of active cores. On these systems // you should specify number of cores manually at compile time. // Without it ALGLIB will run in single-threaded mode. // // Now, back to our example. In this example we will show you: // * how to call SMP version of rmatrixgemm(). Because we work with tiny 2x2 // matrices, we won't expect to see ANY speedup from using multithreading. // The only purpose of this demo is to show how to call SMP functions. // * how to modify number of worker threads used by ALGLIB // double[,] a = new double[,]{{2,1},{1,3}}; double[,] b = new double[,]{{2,1},{0,1}}; double[,] c = new double[,]{{0,0},{0,0}}; int m = 2; int n = 2; int k = 2; double alpha = 1.0; int ia = 0; int ja = 0; int optypea = 0; int ib = 0; int jb = 0; int optypeb = 0; double beta = 0.0; int ic = 0; int jc = 0; // serial code c = new double[,]{{0,0},{0,0}}; alglib.rmatrixgemm(m, n, k, alpha, a, ia, ja, optypea, b, ib, jb, optypeb, beta, ref c, ic, jc); // SMP code with default number of worker threads c = new double[,]{{0,0},{0,0}}; alglib.smp_rmatrixgemm(m, n, k, alpha, a, ia, ja, optypea, b, ib, jb, optypeb, beta, ref c, ic, jc); System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[4,3],[2,4]] // override number of worker threads - use two cores alglib.setnworkers(+2); c = new double[,]{{0,0},{0,0}}; alglib.smp_rmatrixgemm(m, n, k, alpha, a, ia, ja, optypea, b, ib, jb, optypeb, beta, ref c, ic, jc); System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[4,3],[2,4]] System.Console.ReadLine(); return 0; }
public static int Main(string[] args) { // // In this example we assume that you already know how to work with // rmatrixsyrk() function. Below we concentrate on its multithreading // capabilities. // // SMP edition of ALGLIB includes smp_rmatrixsyrk() - multithreaded // version of rmatrixsyrk() function. In the basic edition of ALGLIB // (GPL edition or commercial version without SMP support) this function // just calls single-threaded stub. So, you may call this function from // ANY edition of ALGLIB, but only in SMP edition it will work in really // multithreaded mode. // // In order to use multithreading, you have to: // 1) Install SMP edition of ALGLIB. // 2) This step is specific for C++ users: you should activate OS-specific // capabilities of ALGLIB by defining AE_OS=AE_POSIX (for *nix systems) // or AE_OS=AE_WINDOWS (for Windows systems). // C# users do not have to perform this step because C# programs are // portable across different systems without OS-specific tuning. // 3) Allow ALGLIB to know about number of worker threads to use: // a) autodetection (C++, C#): // ALGLIB will automatically determine number of CPU cores and // (by default) will use all cores except for one. Say, on 4-core // system it will use three cores - unless you manually told it // to use more or less. It will keep your system responsive during // lengthy computations. // Such behavior may be changed with setnworkers() call: // * alglib::setnworkers(0) = use all cores // * alglib::setnworkers(-1) = leave one core unused // * alglib::setnworkers(-2) = leave two cores unused // * alglib::setnworkers(+2) = use 2 cores (even if you have more) // b) manual specification (C++, C#): // You may want to specify maximum number of worker threads during // compile time by means of preprocessor definition AE_NWORKERS. // For C++ it will be "AE_NWORKERS=X" where X can be any positive number. // For C# it is "AE_NWORKERSX", where X should be replaced by number of // workers (AE_NWORKERS2, AE_NWORKERS3, AE_NWORKERS4, ...). // You can add this definition to compiler command line or change // corresponding project settings in your IDE. // // After you installed and configured SMP edition of ALGLIB, you may choose // between serial and multithreaded versions of SMP-capable functions: // * serial version works as usual, in the context of the calling thread // * multithreaded version (with "smp_" prefix) creates (or wakes up) worker // threads, inserts task in the worker queue, and waits for completion of // the task. All processing is done in context of worker thread(s). // // NOTE: because starting/stopping worker threads costs thousands of CPU cycles, // you should not use multithreading for lightweight computational problems. // // NOTE: some old POSIX-compatible operating systems do not support // sysconf(_SC_NPROCESSORS_ONLN) system call which is required in order // to automatically determine number of active cores. On these systems // you should specify number of cores manually at compile time. // Without it ALGLIB will run in single-threaded mode. // // Now, back to our example. In this example we will show you: // * how to call SMP version of rmatrixsyrk(). Because we work with tiny 2x2 // matrices, we won't expect to see ANY speedup from using multithreading. // The only purpose of this demo is to show how to call SMP functions. // * how to modify number of worker threads used by ALGLIB // int n = 2; int k = 1; double alpha = 1.0; int ia = 0; int ja = 0; int optypea = 2; double beta = 0.0; int ic = 0; int jc = 0; bool isupper = true; double[,] a = new double[,]{{1,2}}; double[,] c = new double[0,0]; // // Default number of worker threads. // Preallocate space to store result, call multithreaded version, test. // // NOTE: this function updates only one triangular part of C. In our // example we choose to update upper triangle. // c = new double[,]{{0,0},{0,0}}; alglib.smp_rmatrixsyrk(n, k, alpha, a, ia, ja, optypea, beta, ref c, ic, jc, isupper); System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[1,2],[0,4]] // // Override default number of worker threads (set to 2). // Preallocate space to store result, call multithreaded version, test. // // NOTE: this function updates only one triangular part of C. In our // example we choose to update upper triangle. // alglib.setnworkers(+2); c = new double[,]{{0,0},{0,0}}; alglib.smp_rmatrixsyrk(n, k, alpha, a, ia, ja, optypea, beta, ref c, ic, jc, isupper); System.Console.WriteLine("{0}", alglib.ap.format(c,3)); // EXPECTED: [[1,2],[0,4]] System.Console.ReadLine(); return 0; }
airyf
subpackageairy
function/************************************************************************* Airy function Solution of the differential equation y"(x) = xy. The function returns the two independent solutions Ai, Bi and their first derivatives Ai'(x), Bi'(x). Evaluation is by power series summation for small x, by rational minimax approximations for large x. ACCURACY: Error criterion is absolute when function <= 1, relative when function > 1, except * denotes relative error criterion. For large negative x, the absolute error increases as x^1.5. For large positive x, the relative error increases as x^1.5. Arithmetic domain function # trials peak rms IEEE -10, 0 Ai 10000 1.6e-15 2.7e-16 IEEE 0, 10 Ai 10000 2.3e-14* 1.8e-15* IEEE -10, 0 Ai' 10000 4.6e-15 7.6e-16 IEEE 0, 10 Ai' 10000 1.8e-14* 1.5e-15* IEEE -10, 10 Bi 30000 4.2e-15 5.3e-16 IEEE -10, 10 Bi' 30000 4.9e-15 7.3e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier *************************************************************************/public static void alglib.airy( double x, out double ai, out double aip, out double bi, out double bip)
autogk
subpackage| autogk_d1 | Integrating f=exp(x) by adaptive integrator |
autogkreport
class/************************************************************************* Integration report: * TerminationType = completetion code: * -5 non-convergence of Gauss-Kronrod nodes calculation subroutine. * -1 incorrect parameters were specified * 1 OK * Rep.NFEV countains number of function calculations * Rep.NIntervals contains number of intervals [a,b] was partitioned into. *************************************************************************/public class autogkreport { public int terminationtype; public int nfev; public int nintervals; }
autogkstate
class/************************************************************************* This structure stores state of the integration algorithm. Although this class has public fields, they are not intended for external use. You should use ALGLIB functions to work with this class: * autogksmooth()/AutoGKSmoothW()/... to create objects * autogkintegrate() to begin integration * autogkresults() to get results *************************************************************************/public class autogkstate { }
autogkintegrate
function/************************************************************************* This function is used to launcn iterations of ODE solver It accepts following parameters: diff - callback which calculates dy/dx for given y and x obj - optional object which is passed to diff; can be NULL -- ALGLIB -- Copyright 07.05.2009 by Bochkanov Sergey *************************************************************************/public static void autogkintegrate(autogkstate state, integrator1_func func, object obj)
Examples: [1]
autogkresults
function/************************************************************************* Adaptive integration results Called after AutoGKIteration returned False. Input parameters: State - algorithm state (used by AutoGKIteration). Output parameters: V - integral(f(x)dx,a,b) Rep - optimization report (see AutoGKReport description) -- ALGLIB -- Copyright 14.11.2007 by Bochkanov Sergey *************************************************************************/public static void alglib.autogkresults( autogkstate state, out double v, out autogkreport rep)
Examples: [1]
autogksingular
function/************************************************************************* Integration on a finite interval [A,B]. Integrand have integrable singularities at A/B. F(X) must diverge as "(x-A)^alpha" at A, as "(B-x)^beta" at B, with known alpha/beta (alpha>-1, beta>-1). If alpha/beta are not known, estimates from below can be used (but these estimates should be greater than -1 too). One of alpha/beta variables (or even both alpha/beta) may be equal to 0, which means than function F(x) is non-singular at A/B. Anyway (singular at bounds or not), function F(x) is supposed to be continuous on (A,B). Fast-convergent algorithm based on a Gauss-Kronrod formula is used. Result is calculated with accuracy close to the machine precision. INPUT PARAMETERS: A, B - interval boundaries (A<B, A=B or A>B) Alpha - power-law coefficient of the F(x) at A, Alpha>-1 Beta - power-law coefficient of the F(x) at B, Beta>-1 OUTPUT PARAMETERS State - structure which stores algorithm state SEE ALSO AutoGKSmooth, AutoGKSmoothW, AutoGKResults. -- ALGLIB -- Copyright 06.05.2009 by Bochkanov Sergey *************************************************************************/public static void alglib.autogksingular( double a, double b, double alpha, double beta, out autogkstate state)
autogksmooth
function/************************************************************************* Integration of a smooth function F(x) on a finite interval [a,b]. Fast-convergent algorithm based on a Gauss-Kronrod formula is used. Result is calculated with accuracy close to the machine precision. Algorithm works well only with smooth integrands. It may be used with continuous non-smooth integrands, but with less performance. It should never be used with integrands which have integrable singularities at lower or upper limits - algorithm may crash. Use AutoGKSingular in such cases. INPUT PARAMETERS: A, B - interval boundaries (A<B, A=B or A>B) OUTPUT PARAMETERS State - structure which stores algorithm state SEE ALSO AutoGKSmoothW, AutoGKSingular, AutoGKResults. -- ALGLIB -- Copyright 06.05.2009 by Bochkanov Sergey *************************************************************************/public static void alglib.autogksmooth( double a, double b, out autogkstate state)
Examples: [1]
autogksmoothw
function/************************************************************************* Integration of a smooth function F(x) on a finite interval [a,b]. This subroutine is same as AutoGKSmooth(), but it guarantees that interval [a,b] is partitioned into subintervals which have width at most XWidth. Subroutine can be used when integrating nearly-constant function with narrow "bumps" (about XWidth wide). If "bumps" are too narrow, AutoGKSmooth subroutine can overlook them. INPUT PARAMETERS: A, B - interval boundaries (A<B, A=B or A>B) OUTPUT PARAMETERS State - structure which stores algorithm state SEE ALSO AutoGKSmooth, AutoGKSingular, AutoGKResults. -- ALGLIB -- Copyright 06.05.2009 by Bochkanov Sergey *************************************************************************/public static void alglib.autogksmoothw( double a, double b, double xwidth, out autogkstate state)
public static void int_function_1_func(double x, double xminusa, double bminusx, ref double y, object obj) { // this callback calculates f(x)=exp(x) y = Math.Exp(x); } public static int Main(string[] args) { // // This example demonstrates integration of f=exp(x) on [0,1]: // * first, autogkstate is initialized // * then we call integration function // * and finally we obtain results with autogkresults() call // double a = 0; double b = 1; alglib.autogkstate s; double v; alglib.autogkreport rep; alglib.autogksmooth(a, b, out s); alglib.autogkintegrate(s, int_function_1_func, null); alglib.autogkresults(s, out v, out rep); System.Console.WriteLine("{0:F2", v); // EXPECTED: 1.7182 System.Console.ReadLine(); return 0; }
basestat
subpackage| basestat_d_base | Basic functionality (moments, adev, median, percentile) | |
| basestat_d_c2 | Correlation (covariance) between two random variables | |
| basestat_d_cm | Correlation (covariance) between components of random vector | |
| basestat_d_cm2 | Correlation (covariance) between two random vectors |
cov2
function/************************************************************************* 2-sample covariance Input parameters: X - sample 1 (array indexes: [0..N-1]) Y - sample 2 (array indexes: [0..N-1]) N - N>=0, sample size: * if given, only N leading elements of X/Y are processed * if not given, automatically determined from input sizes Result: covariance (zero for N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static double alglib.cov2(double[] x, double[] y) public static double alglib.cov2(double[] x, double[] y, int n)
Examples: [1]
covm
function/************************************************************************* Covariance matrix SMP EDITION OF ALGLIB: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Although multicore version is pretty efficient on ! large problems, we do not recommend you to use it on small problems - ! with covariance matrices smaller than 128*128. INPUT PARAMETERS: X - array[N,M], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation N - N>=0, number of observations: * if given, only leading N rows of X are used * if not given, automatically determined from input size M - M>0, number of variables: * if given, only leading M columns of X are used * if not given, automatically determined from input size OUTPUT PARAMETERS: C - array[M,M], covariance matrix (zero if N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static void alglib.covm(double[,] x, out double[,] c) public static void alglib.covm(double[,] x, int n, int m, out double[,] c) public static void alglib.smp_covm(double[,] x, out double[,] c) public static void alglib.smp_covm(double[,] x, int n, int m, out double[,] c)
Examples: [1]
covm2
function/************************************************************************* Cross-covariance matrix SMP EDITION OF ALGLIB: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Although multicore version is pretty efficient on ! large problems, we do not recommend you to use it on small problems - ! with covariance matrices smaller than 128*128. INPUT PARAMETERS: X - array[N,M1], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation Y - array[N,M2], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation N - N>=0, number of observations: * if given, only leading N rows of X/Y are used * if not given, automatically determined from input sizes M1 - M1>0, number of variables in X: * if given, only leading M1 columns of X are used * if not given, automatically determined from input size M2 - M2>0, number of variables in Y: * if given, only leading M1 columns of X are used * if not given, automatically determined from input size OUTPUT PARAMETERS: C - array[M1,M2], cross-covariance matrix (zero if N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static void alglib.covm2(double[,] x, double[,] y, out double[,] c) public static void alglib.covm2( double[,] x, double[,] y, int n, int m1, int m2, out double[,] c) public static void alglib.smp_covm2(double[,] x, double[,] y, out double[,] c) public static void alglib.smp_covm2( double[,] x, double[,] y, int n, int m1, int m2, out double[,] c)
Examples: [1]
pearsoncorr2
function/************************************************************************* Pearson product-moment correlation coefficient Input parameters: X - sample 1 (array indexes: [0..N-1]) Y - sample 2 (array indexes: [0..N-1]) N - N>=0, sample size: * if given, only N leading elements of X/Y are processed * if not given, automatically determined from input sizes Result: Pearson product-moment correlation coefficient (zero for N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static double alglib.pearsoncorr2(double[] x, double[] y) public static double alglib.pearsoncorr2(double[] x, double[] y, int n)
Examples: [1]
pearsoncorrelation
function/************************************************************************* Obsolete function, we recommend to use PearsonCorr2(). -- ALGLIB -- Copyright 09.04.2007 by Bochkanov Sergey *************************************************************************/public static double alglib.pearsoncorrelation( double[] x, double[] y, int n)
pearsoncorrm
function/************************************************************************* Pearson product-moment correlation matrix SMP EDITION OF ALGLIB: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Although multicore version is pretty efficient on ! large problems, we do not recommend you to use it on small problems - ! with correlation matrices smaller than 128*128. INPUT PARAMETERS: X - array[N,M], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation N - N>=0, number of observations: * if given, only leading N rows of X are used * if not given, automatically determined from input size M - M>0, number of variables: * if given, only leading M columns of X are used * if not given, automatically determined from input size OUTPUT PARAMETERS: C - array[M,M], correlation matrix (zero if N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static void alglib.pearsoncorrm(double[,] x, out double[,] c) public static void alglib.pearsoncorrm( double[,] x, int n, int m, out double[,] c) public static void alglib.smp_pearsoncorrm(double[,] x, out double[,] c) public static void alglib.smp_pearsoncorrm( double[,] x, int n, int m, out double[,] c)
Examples: [1]
pearsoncorrm2
function/************************************************************************* Pearson product-moment cross-correlation matrix SMP EDITION OF ALGLIB: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Although multicore version is pretty efficient on ! large problems, we do not recommend you to use it on small problems - ! with correlation matrices smaller than 128*128. INPUT PARAMETERS: X - array[N,M1], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation Y - array[N,M2], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation N - N>=0, number of observations: * if given, only leading N rows of X/Y are used * if not given, automatically determined from input sizes M1 - M1>0, number of variables in X: * if given, only leading M1 columns of X are used * if not given, automatically determined from input size M2 - M2>0, number of variables in Y: * if given, only leading M1 columns of X are used * if not given, automatically determined from input size OUTPUT PARAMETERS: C - array[M1,M2], cross-correlation matrix (zero if N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static void alglib.pearsoncorrm2( double[,] x, double[,] y, out double[,] c) public static void alglib.pearsoncorrm2( double[,] x, double[,] y, int n, int m1, int m2, out double[,] c) public static void alglib.smp_pearsoncorrm2( double[,] x, double[,] y, out double[,] c) public static void alglib.smp_pearsoncorrm2( double[,] x, double[,] y, int n, int m1, int m2, out double[,] c)
Examples: [1]
rankdata
function/************************************************************************* *************************************************************************/public static void alglib.rankdata(ref double[,] xy) public static void alglib.rankdata( ref double[,] xy, int npoints, int nfeatures) public static void alglib.smp_rankdata(ref double[,] xy) public static void alglib.smp_rankdata( ref double[,] xy, int npoints, int nfeatures)
rankdatacentered
function/************************************************************************* *************************************************************************/public static void alglib.rankdatacentered(ref double[,] xy) public static void alglib.rankdatacentered( ref double[,] xy, int npoints, int nfeatures) public static void alglib.smp_rankdatacentered(ref double[,] xy) public static void alglib.smp_rankdatacentered( ref double[,] xy, int npoints, int nfeatures)
sampleadev
function/************************************************************************* ADev Input parameters: X - sample N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X Output parameters: ADev- ADev -- ALGLIB -- Copyright 06.09.2006 by Bochkanov Sergey *************************************************************************/public static void alglib.sampleadev(double[] x, out double adev) public static void alglib.sampleadev(double[] x, int n, out double adev)
Examples: [1]
samplekurtosis
function/************************************************************************* Calculation of the kurtosis. INPUT PARAMETERS: X - sample N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X NOTE: This function return result which calculated by 'SampleMoments' function and stored at 'Kurtosis' variable. -- ALGLIB -- Copyright 06.09.2006 by Bochkanov Sergey *************************************************************************/public static double alglib.samplekurtosis(double[] x) public static double alglib.samplekurtosis(double[] x, int n)
samplemean
function/************************************************************************* Calculation of the mean. INPUT PARAMETERS: X - sample N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X NOTE: This function return result which calculated by 'SampleMoments' function and stored at 'Mean' variable. -- ALGLIB -- Copyright 06.09.2006 by Bochkanov Sergey *************************************************************************/public static double alglib.samplemean(double[] x) public static double alglib.samplemean(double[] x, int n)
samplemedian
function/************************************************************************* Median calculation. Input parameters: X - sample (array indexes: [0..N-1]) N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X Output parameters: Median -- ALGLIB -- Copyright 06.09.2006 by Bochkanov Sergey *************************************************************************/public static void alglib.samplemedian(double[] x, out double median) public static void alglib.samplemedian( double[] x, int n, out double median)
Examples: [1]
samplemoments
function/************************************************************************* Calculation of the distribution moments: mean, variance, skewness, kurtosis. INPUT PARAMETERS: X - sample N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X OUTPUT PARAMETERS Mean - mean. Variance- variance. Skewness- skewness (if variance<>0; zero otherwise). Kurtosis- kurtosis (if variance<>0; zero otherwise).  Input parameters: X - sample (array indexes: [0..N-1]) N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X P - percentile (0<=P<=1) Output parameters: V - percentile -- ALGLIB -- Copyright 01.03.2008 by Bochkanov Sergey *************************************************************************/public static void alglib.samplepercentile( double[] x, double p, out double v) public static void alglib.samplepercentile( double[] x, int n, double p, out double v)
Examples: [1]
sampleskewness
function/************************************************************************* Calculation of the skewness. INPUT PARAMETERS: X - sample N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X NOTE: This function return result which calculated by 'SampleMoments' function and stored at 'Skewness' variable. -- ALGLIB -- Copyright 06.09.2006 by Bochkanov Sergey *************************************************************************/public static double alglib.sampleskewness(double[] x) public static double alglib.sampleskewness(double[] x, int n)
samplevariance
function/************************************************************************* Calculation of the variance. INPUT PARAMETERS: X - sample N - N>=0, sample size: * if given, only leading N elements of X are processed * if not given, automatically determined from size of X NOTE: This function return result which calculated by 'SampleMoments' function and stored at 'Variance' variable. -- ALGLIB -- Copyright 06.09.2006 by Bochkanov Sergey *************************************************************************/public static double alglib.samplevariance(double[] x) public static double alglib.samplevariance(double[] x, int n)
spearmancorr2
function/************************************************************************* Spearman's rank correlation coefficient Input parameters: X - sample 1 (array indexes: [0..N-1]) Y - sample 2 (array indexes: [0..N-1]) N - N>=0, sample size: * if given, only N leading elements of X/Y are processed * if not given, automatically determined from input sizes Result: Spearman's rank correlation coefficient (zero for N=0 or N=1) -- ALGLIB -- Copyright 09.04.2007 by Bochkanov Sergey *************************************************************************/public static double alglib.spearmancorr2(double[] x, double[] y) public static double alglib.spearmancorr2(double[] x, double[] y, int n)
Examples: [1]
spearmancorrm
function/************************************************************************* Spearman's rank correlation matrix SMP EDITION OF ALGLIB: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Although multicore version is pretty efficient on ! large problems, we do not recommend you to use it on small problems - ! with correlation matrices smaller than 128*128. INPUT PARAMETERS: X - array[N,M], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation N - N>=0, number of observations: * if given, only leading N rows of X are used * if not given, automatically determined from input size M - M>0, number of variables: * if given, only leading M columns of X are used * if not given, automatically determined from input size OUTPUT PARAMETERS: C - array[M,M], correlation matrix (zero if N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static void alglib.spearmancorrm(double[,] x, out double[,] c) public static void alglib.spearmancorrm( double[,] x, int n, int m, out double[,] c) public static void alglib.smp_spearmancorrm(double[,] x, out double[,] c) public static void alglib.smp_spearmancorrm( double[,] x, int n, int m, out double[,] c)
Examples: [1]
spearmancorrm2
function/************************************************************************* Spearman's rank cross-correlation matrix SMP EDITION OF ALGLIB: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Although multicore version is pretty efficient on ! large problems, we do not recommend you to use it on small problems - ! with correlation matrices smaller than 128*128. INPUT PARAMETERS: X - array[N,M1], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation Y - array[N,M2], sample matrix: * J-th column corresponds to J-th variable * I-th row corresponds to I-th observation N - N>=0, number of observations: * if given, only leading N rows of X/Y are used * if not given, automatically determined from input sizes M1 - M1>0, number of variables in X: * if given, only leading M1 columns of X are used * if not given, automatically determined from input size M2 - M2>0, number of variables in Y: * if given, only leading M1 columns of X are used * if not given, automatically determined from input size OUTPUT PARAMETERS: C - array[M1,M2], cross-correlation matrix (zero if N=0 or N=1) -- ALGLIB -- Copyright 28.10.2010 by Bochkanov Sergey *************************************************************************/public static void alglib.spearmancorrm2( double[,] x, double[,] y, out double[,] c) public static void alglib.spearmancorrm2( double[,] x, double[,] y, int n, int m1, int m2, out double[,] c) public static void alglib.smp_spearmancorrm2( double[,] x, double[,] y, out double[,] c) public static void alglib.smp_spearmancorrm2( double[,] x, double[,] y, int n, int m1, int m2, out double[,] c)
Examples: [1]
spearmanrankcorrelation
function/************************************************************************* Obsolete function, we recommend to use SpearmanCorr2(). -- ALGLIB -- Copyright 09.04.2007 by Bochkanov Sergey *************************************************************************/public static double alglib.spearmanrankcorrelation( double[] x, double[] y, int n)
public static int Main(string[] args) { double[] x = new double[]{0,1,4,9,16,25,36,49,64,81}; double mean; double variance; double skewness; double kurtosis; double adev; double p; double v; // // Here we demonstrate calculation of sample moments // (mean, variance, skewness, kurtosis) // alglib.samplemoments(x, out mean, out variance, out skewness, out kurtosis); System.Console.WriteLine("{0:F1", mean); // EXPECTED: 28.5 System.Console.WriteLine("{0:F1", variance); // EXPECTED: 801.1667 System.Console.WriteLine("{0:F1", skewness); // EXPECTED: 0.5751 System.Console.WriteLine("{0:F1", kurtosis); // EXPECTED: -1.2666 // // Average deviation // alglib.sampleadev(x, out adev); System.Console.WriteLine("{0:F1", adev); // EXPECTED: 23.2 // // Median and percentile // alglib.samplemedian(x, out v); System.Console.WriteLine("{0:F1", v); // EXPECTED: 20.5 p = 0.5; alglib.samplepercentile(x, p, out v); System.Console.WriteLine("{0:F1", v); // EXPECTED: 20.5 System.Console.ReadLine(); return 0; }
public static int Main(string[] args) { // // We have two samples - x and y, and want to measure dependency between them // double[] x = new double[]{0,1,4,9,16,25,36,49,64,81}; double[] y = new double[]{0,1,2,3,4,5,6,7,8,9}; double v; // // Three dependency measures are calculated: // * covariation // * Pearson correlation // * Spearman rank correlation // v = alglib.cov2(x, y); System.Console.WriteLine("{0:F2", v); // EXPECTED: 82.5 v = alglib.pearsoncorr2(x, y); System.Console.WriteLine("{0:F2", v); // EXPECTED: 0.9627 v = alglib.spearmancorr2(x, y); System.Console.WriteLine("{0:F2", v); // EXPECTED: 1.000 System.Console.ReadLine(); return 0; }
public static int Main(string[] args) { // // X is a sample matrix: // * I-th row corresponds to I-th observation // * J-th column corresponds to J-th variable // double[,] x = new double[,]{{1,0,1},{1,1,0},{-1,1,0},{-2,-1,1},{-1,0,9}}; double[,] c; // // Three dependency measures are calculated: // * covariation // * Pearson correlation // * Spearman rank correlation // // Result is stored into C, with C[i,j] equal to correlation // (covariance) between I-th and J-th variables of X. // alglib.covm(x, out c); System.Console.WriteLine("{0}", alglib.ap.format(c,1)); // EXPECTED: [[1.80,0.60,-1.40],[0.60,0.70,-0.80],[-1.40,-0.80,14.70]] alglib.pearsoncorrm(x, out c); System.Console.WriteLine("{0}", alglib.ap.format(c,1)); // EXPECTED: [[1.000,0.535,-0.272],[0.535,1.000,-0.249],[-0.272,-0.249,1.000]] alglib.spearmancorrm(x, out c); System.Console.WriteLine("{0}", alglib.ap.format(c,1)); // EXPECTED: [[1.000,0.556,-0.306],[0.556,1.000,-0.750],[-0.306,-0.750,1.000]] System.Console.ReadLine(); return 0; }
public static int Main(string[] args) { // // X and Y are sample matrices: // * I-th row corresponds to I-th observation // * J-th column corresponds to J-th variable // double[,] x = new double[,]{{1,0,1},{1,1,0},{-1,1,0},{-2,-1,1},{-1,0,9}}; double[,] y = new double[,]{{2,3},{2,1},{-1,6},{-9,9},{7,1}}; double[,] c; // // Three dependency measures are calculated: // * covariation // * Pearson correlation // * Spearman rank correlation // // Result is stored into C, with C[i,j] equal to correlation // (covariance) between I-th variable of X and J-th variable of Y. // alglib.covm2(x, y, out c); System.Console.WriteLine("{0}", alglib.ap.format(c,1)); // EXPECTED: [[4.100,-3.250],[2.450,-1.500],[13.450,-5.750]] alglib.pearsoncorrm2(x, y, out c); System.Console.WriteLine("{0}", alglib.ap.format(c,1)); // EXPECTED: [[0.519,-0.699],[0.497,-0.518],[0.596,-0.433]] alglib.spearmancorrm2(x, y, out c); System.Console.WriteLine("{0}", alglib.ap.format(c,1)); // EXPECTED: [[0.541,-0.649],[0.216,-0.433],[0.433,-0.135]] System.Console.ReadLine(); return 0; }
bdss
subpackagedsoptimalsplit2
function/************************************************************************* Optimal binary classification Algorithms finds optimal (=with minimal cross-entropy) binary partition. Internal subroutine. INPUT PARAMETERS: A - array[0..N-1], variable C - array[0..N-1], class numbers (0 or 1). N - array size OUTPUT PARAMETERS: Info - completetion code: * -3, all values of A[] are same (partition is impossible) * -2, one of C[] is incorrect (<0, >1) * -1, incorrect pararemets were passed (N<=0). * 1, OK Threshold- partiton boundary. Left part contains values which are strictly less than Threshold. Right part contains values which are greater than or equal to Threshold. PAL, PBL- probabilities P(0|v<Threshold) and P(1|v<Threshold) PAR, PBR- probabilities P(0|v>=Threshold) and P(1|v>=Threshold) CVE - cross-validation estimate of cross-entropy -- ALGLIB -- Copyright 22.05.2008 by Bochkanov Sergey *************************************************************************/public static void alglib.dsoptimalsplit2( double[] a, int[] c, int n, out int info, out double threshold, out double pal, out double pbl, out double par, out double pbr, out double cve)
dsoptimalsplit2fast
function/************************************************************************* Optimal partition, internal subroutine. Fast version. Accepts: A array[0..N-1] array of attributes array[0..N-1] C array[0..N-1] array of class labels TiesBuf array[0..N] temporaries (ties) CntBuf array[0..2*NC-1] temporaries (counts) Alpha centering factor (0<=alpha<=1, recommended value - 0.05) BufR array[0..N-1] temporaries BufI array[0..N-1] temporaries Output: Info error code (">0"=OK, "<0"=bad) RMS training set RMS error CVRMS leave-one-out RMS error Note: content of all arrays is changed by subroutine; it doesn't allocate temporaries. -- ALGLIB -- Copyright 11.12.2008 by Bochkanov Sergey *************************************************************************/public static void alglib.dsoptimalsplit2fast( ref double[] a, ref int[] c, ref int[] tiesbuf, ref int[] cntbuf, ref double[] bufr, ref int[] bufi, int n, int nc, double alpha, out int info, out double threshold, out double rms, out double cvrms)
bdsvd
subpackagermatrixbdsvd
function/************************************************************************* Singular value decomposition of a bidiagonal matrix (extended algorithm) The algorithm performs the singular value decomposition of a bidiagonal matrix B (upper or lower) representing it as B = Q*S*P^T, where Q and P - orthogonal matrices, S - diagonal matrix with non-negative elements on the main diagonal, in descending order. The algorithm finds singular values. In addition, the algorithm can calculate matrices Q and P (more precisely, not the matrices, but their product with given matrices U and VT - U*Q and (P^T)*VT)). Of course, matrices U and VT can be of any type, including identity. Furthermore, the algorithm can calculate Q'*C (this product is calculated more effectively than U*Q, because this calculation operates with rows instead of matrix columns). The feature of the algorithm is its ability to find all singular values including those which are arbitrarily close to 0 with relative accuracy close to machine precision. If the parameter IsFractionalAccuracyRequired is set to True, all singular values will have high relative accuracy close to machine precision. If the parameter is set to False, only the biggest singular value will have relative accuracy close to machine precision. The absolute error of other singular values is equal to the absolute error of the biggest singular value. Input parameters: D - main diagonal of matrix B. Array whose index ranges within [0..N-1]. E - superdiagonal (or subdiagonal) of matrix B. Array whose index ranges within [0..N-2]. N - size of matrix B. IsUpper - True, if the matrix is upper bidiagonal.  U - matrix to be multiplied by Q. Array whose indexes range within [0..NRU-1, 0..N-1]. The matrix can be bigger, in that case only the submatrix [0..NRU-1, 0..N-1] will be multiplied by Q. NRU - number of rows in matrix U. C - matrix to be multiplied by Q'. Array whose indexes range within [0..N-1, 0..NCC-1]. The matrix can be bigger, in that case only the submatrix [0..N-1, 0..NCC-1] will be multiplied by Q'. NCC - number of columns in matrix C. VT - matrix to be multiplied by P^T. Array whose indexes range within [0..N-1, 0..NCVT-1]. The matrix can be bigger, in that case only the submatrix [0..N-1, 0..NCVT-1] will be multiplied by P^T. NCVT - number of columns in matrix VT. Output parameters: D - singular values of matrix B in descending order. U - if NRU>0, contains matrix U*Q. VT - if NCVT>0, contains matrix (P^T)*VT. C - if NCC>0, contains matrix Q'*C. Result: True, if the algorithm has converged. False, if the algorithm hasn't converged (rare case). Additional information: The type of convergence is controlled by the internal parameter TOL. If the parameter is greater than 0, the singular values will have relative accuracy TOL. If TOL<0, the singular values will have absolute accuracy ABS(TOL)*norm(B). By default, |TOL| falls within the range of 10*Epsilon and 100*Epsilon, where Epsilon is the machine precision. It is not recommended to use TOL less than 10*Epsilon since this will considerably slow down the algorithm and may not lead to error decreasing. History: * 31 March, 2007. changed MAXITR from 6 to 12. -- LAPACK routine (version 3.0) -- Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., Courant Institute, Argonne National Lab, and Rice University October 31, 1999. *************************************************************************/public static bool alglib.rmatrixbdsvd( ref double[] d, double[] e, int n, bool isupper, bool isfractionalaccuracyrequired, ref double[,] u, int nru, ref double[,] c, int ncc, ref double[,] vt, int ncvt)
bessel
subpackagebesseli0
function/************************************************************************* Modified Bessel function of order zero Returns modified Bessel function of order zero of the argument. The function is defined as i0(x) = j0( ix ). The range is partitioned into the two intervals [0,8] and (8, infinity). Chebyshev polynomial expansions are employed in each interval. ACCURACY: Relative error: arithmetic domain # trials peak rms IEEE 0,30 30000 5.8e-16 1.4e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besseli0(double x)
besseli1
function/************************************************************************* Modified Bessel function of order one Returns modified Bessel function of order one of the argument. The function is defined as i1(x) = -i j1( ix ). The range is partitioned into the two intervals [0,8] and (8, infinity). Chebyshev polynomial expansions are employed in each interval. ACCURACY: Relative error: arithmetic domain # trials peak rms IEEE 0, 30 30000 1.9e-15 2.1e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1985, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besseli1(double x)
besselj0
function/************************************************************************* Bessel function of order zero Returns Bessel function of order zero of the argument. The domain is divided into the intervals [0, 5] and (5, infinity). In the first interval the following rational approximation is used: 2 2 (w - r ) (w - r ) P (w) / Q (w) 1 2 3 8 2 where w = x and the two r's are zeros of the function. In the second interval, the Hankel asymptotic expansion is employed with two rational functions of degree 6/6 and 7/7. ACCURACY: Absolute error: arithmetic domain # trials peak rms IEEE 0, 30 60000 4.2e-16 1.1e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besselj0(double x)
besselj1
function/************************************************************************* Bessel function of order one Returns Bessel function of order one of the argument. The domain is divided into the intervals [0, 8] and (8, infinity). In the first interval a 24 term Chebyshev expansion is used. In the second, the asymptotic trigonometric representation is employed using two rational functions of degree 5/5. ACCURACY: Absolute error: arithmetic domain # trials peak rms IEEE 0, 30 30000 2.6e-16 1.1e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besselj1(double x)
besseljn
function/************************************************************************* Bessel function of integer order Returns Bessel function of order n, where n is a (possibly negative) integer. The ratio of jn(x) to j0(x) is computed by backward recurrence. First the ratio jn/jn-1 is found by a continued fraction expansion. Then the recurrence relating successive orders is applied until j0 or j1 is reached. If n = 0 or 1 the routine for j0 or j1 is called directly. ACCURACY: Absolute error: arithmetic range # trials peak rms IEEE 0, 30 5000 4.4e-16 7.9e-17 Not suitable for large n or x. Use jv() (fractional order) instead. Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besseljn(int n, double x)
besselk0
function/************************************************************************* Modified Bessel function, second kind, order zero Returns modified Bessel function of the second kind of order zero of the argument. The range is partitioned into the two intervals [0,8] and (8, infinity). Chebyshev polynomial expansions are employed in each interval. ACCURACY: Tested at 2000 random points between 0 and 8. Peak absolute error (relative when K0 > 1) was 1.46e-14; rms, 4.26e-15. Relative error: arithmetic domain # trials peak rms IEEE 0, 30 30000 1.2e-15 1.6e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besselk0(double x)
besselk1
function/************************************************************************* Modified Bessel function, second kind, order one Computes the modified Bessel function of the second kind of order one of the argument. The range is partitioned into the two intervals [0,2] and (2, infinity). Chebyshev polynomial expansions are employed in each interval. ACCURACY: Relative error: arithmetic domain # trials peak rms IEEE 0, 30 30000 1.2e-15 1.6e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besselk1(double x)
besselkn
function/************************************************************************* Modified Bessel function, second kind, integer order Returns modified Bessel function of the second kind of order n of the argument. The range is partitioned into the two intervals [0,9.55] and (9.55, infinity). An ascending power series is used in the low range, and an asymptotic expansion in the high range. ACCURACY: Relative error: arithmetic domain # trials peak rms IEEE 0,30 90000 1.8e-8 3.0e-10 Error is high only near the crossover point x = 9.55 between the two expansions used. Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1988, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besselkn(int nn, double x)
bessely0
function/************************************************************************* Bessel function of the second kind, order zero Returns Bessel function of the second kind, of order zero, of the argument. The domain is divided into the intervals [0, 5] and (5, infinity). In the first interval a rational approximation R(x) is employed to compute y0(x) = R(x) + 2 * log(x) * j0(x) / PI. Thus a call to j0() is required. In the second interval, the Hankel asymptotic expansion is employed with two rational functions of degree 6/6 and 7/7. ACCURACY: Absolute error, when y0(x) < 1; else relative error: arithmetic domain # trials peak rms IEEE 0, 30 30000 1.3e-15 1.6e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.bessely0(double x)
bessely1
function/************************************************************************* Bessel function of second kind of order one Returns Bessel function of the second kind of order one of the argument. The domain is divided into the intervals [0, 8] and (8, infinity). In the first interval a 25 term Chebyshev expansion is used, and a call to j1() is required. In the second, the asymptotic trigonometric representation is employed using two rational functions of degree 5/5. ACCURACY: Absolute error: arithmetic domain # trials peak rms IEEE 0, 30 30000 1.0e-15 1.3e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1989, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.bessely1(double x)
besselyn
function/************************************************************************* Bessel function of second kind of integer order Returns Bessel function of order n, where n is a (possibly negative) integer. The function is evaluated by forward recurrence on n, starting with values computed by the routines y0() and y1(). If n = 0 or 1 the routine for y0 or y1 is called directly. ACCURACY: Absolute error, except relative when y > 1: arithmetic domain # trials peak rms IEEE 0, 30 30000 3.4e-15 4.3e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.besselyn(int n, double x)
betaf
subpackagebeta
function/************************************************************************* Beta function - - | (a) | (b) beta( a, b ) = -----------. - | (a+b) For large arguments the logarithm of the function is evaluated using lgam(), then exponentiated. ACCURACY: Relative error: arithmetic domain # trials peak rms IEEE 0,30 30000 8.1e-14 1.1e-14 Cephes Math Library Release 2.0: April, 1987 Copyright 1984, 1987 by Stephen L. Moshier *************************************************************************/public static double alglib.beta(double a, double b)
binomialdistr
subpackagebinomialcdistribution
function/************************************************************************* Complemented binomial distribution Returns the sum of the terms k+1 through n of the Binomial probability density: n -- ( n ) j n-j > ( ) p (1-p) -- ( j ) j=k+1 The terms are not summed directly; instead the incomplete beta integral is employed, according to the formula y = bdtrc( k, n, p ) = incbet( k+1, n-k, p ). The arguments must be positive, with p ranging from 0 to 1. ACCURACY: Tested at random points (a,b,p). a,b Relative error: arithmetic domain # trials peak rms For p between 0.001 and 1: IEEE 0,100 100000 6.7e-15 8.2e-16 For p between 0 and .001: IEEE 0,100 100000 1.5e-13 2.7e-15 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1995, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.binomialcdistribution(int k, int n, double p)
binomialdistribution
function/************************************************************************* Binomial distribution Returns the sum of the terms 0 through k of the Binomial probability density: k -- ( n ) j n-j > ( ) p (1-p) -- ( j ) j=0 The terms are not summed directly; instead the incomplete beta integral is employed, according to the formula y = bdtr( k, n, p ) = incbet( n-k, k+1, 1-p ). The arguments must be positive, with p ranging from 0 to 1. ACCURACY: Tested at random points (a,b,p), with p between 0 and 1. a,b Relative error: arithmetic domain # trials peak rms For p between 0.001 and 1: IEEE 0,100 100000 4.3e-15 2.6e-16 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1995, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.binomialdistribution(int k, int n, double p)
invbinomialdistribution
function/************************************************************************* Inverse binomial distribution Finds the event probability p such that the sum of the terms 0 through k of the Binomial probability density is equal to the given cumulative probability y. This is accomplished using the inverse beta integral function and the relation 1 - p = incbi( n-k, k+1, y ). ACCURACY: Tested at random points (a,b,p). a,b Relative error: arithmetic domain # trials peak rms For p between 0.001 and 1: IEEE 0,100 100000 2.3e-14 6.4e-16 IEEE 0,10000 100000 6.6e-12 1.2e-13 For p between 10^-6 and 0.001: IEEE 0,100 100000 2.0e-12 1.3e-14 IEEE 0,10000 100000 1.5e-12 3.2e-14 Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 1995, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.invbinomialdistribution( int k, int n, double y)
chebyshev
subpackagechebyshevcalculate
function/************************************************************************* Calculation of the value of the Chebyshev polynomials of the first and second kinds. Parameters: r - polynomial kind, either 1 or 2. n - degree, n>=0 x - argument, -1 <= x <= 1 Result: the value of the Chebyshev polynomial at x *************************************************************************/public static double alglib.chebyshevcalculate(int r, int n, double x)
chebyshevcoefficients
function/************************************************************************* Representation of Tn as C[0] + C[1]*X + ... + C[N]*X^N Input parameters: N - polynomial degree, n>=0 Output parameters: C - coefficients *************************************************************************/public static void alglib.chebyshevcoefficients(int n, out double[] c)
chebyshevsum
function/************************************************************************* Summation of Chebyshev polynomials using Clenshaw’s recurrence formula. This routine calculates c[0]*T0(x) + c[1]*T1(x) + ... + c[N]*TN(x) or c[0]*U0(x) + c[1]*U1(x) + ... + c[N]*UN(x) depending on the R. Parameters: r - polynomial kind, either 1 or 2. n - degree, n>=0 x - argument Result: the value of the Chebyshev polynomial at x *************************************************************************/public static double alglib.chebyshevsum( double[] c, int r, int n, double x)
fromchebyshev
function/************************************************************************* Conversion of a series of Chebyshev polynomials to a power series. Represents A[0]*T0(x) + A[1]*T1(x) + ... + A[N]*Tn(x) as B[0] + B[1]*X + ... + B[N]*X^N. Input parameters: A - Chebyshev series coefficients N - degree, N>=0 Output parameters B - power series coefficients *************************************************************************/public static void alglib.fromchebyshev(double[] a, int n, out double[] b)
chisquaredistr
subpackagechisquarecdistribution
function/************************************************************************* Complemented Chi-square distribution Returns the area under the right hand tail (from x to infinity) of the Chi square probability density function with v degrees of freedom: inf. - 1 | | v/2-1 -t/2 P( x | v ) = ----------- | t e dt v/2 - | | 2 | (v/2) - x where x is the Chi-square variable. The incomplete gamma integral is used, according to the formula y = chdtr( v, x ) = igamc( v/2.0, x/2.0 ). The arguments must both be positive. ACCURACY: See incomplete gamma function Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.chisquarecdistribution(double v, double x)
chisquaredistribution
function/************************************************************************* Chi-square distribution Returns the area under the left hand tail (from 0 to x) of the Chi square probability density function with v degrees of freedom. x - 1 | | v/2-1 -t/2 P( x | v ) = ----------- | t e dt v/2 - | | 2 | (v/2) - 0 where x is the Chi-square variable. The incomplete gamma integral is used, according to the formula y = chdtr( v, x ) = igam( v/2.0, x/2.0 ). The arguments must both be positive. ACCURACY: See incomplete gamma function Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.chisquaredistribution(double v, double x)
invchisquaredistribution
function/************************************************************************* Inverse of complemented Chi-square distribution Finds the Chi-square argument x such that the integral from x to infinity of the Chi-square density is equal to the given cumulative probability y. This is accomplished using the inverse gamma integral function and the relation x/2 = igami( df/2, y ); ACCURACY: See inverse incomplete gamma function Cephes Math Library Release 2.8: June, 2000 Copyright 1984, 1987, 2000 by Stephen L. Moshier *************************************************************************/public static double alglib.invchisquaredistribution(double v, double y)
clustering
subpackage| clst_ahc | Simple hierarchical clusterization with Euclidean distance function | |
| clst_distance | Clusterization with different metric types | |
| clst_kclusters | Obtaining K top clusters from clusterization tree | |
| clst_kmeans | Simple k-means clusterization | |
| clst_linkage | Clusterization with different linkage types |
ahcreport
class/************************************************************************* This structure is used to store results of the agglomerative hierarchical clustering (AHC). Following information is returned: * NPoints contains number of points in the original dataset * Z contains information about merges performed (see below). Z contains indexes from the original (unsorted) dataset and it can be used when you need to know what points were merged. However, it is not convenient when you want to build a dendrograd (see below). * if you want to build dendrogram, you can use Z, but it is not good option, because Z contains indexes from unsorted dataset. Dendrogram built from such dataset is likely to have intersections. So, you have to reorder you points before building dendrogram. Permutation which reorders point is returned in P. Another representation of merges, which is more convenient for dendorgram construction, is returned in PM. * more information on format of Z, P and PM can be found below and in the examples from ALGLIB Reference Manual. FORMAL DESCRIPTION OF FIELDS: NPoints number of points Z array[NPoints-1,2], contains indexes of clusters linked in pairs to form clustering tree. I-th row corresponds to I-th merge: * Z[I,0] - index of the first cluster to merge * Z[I,1] - index of the second cluster to merge * Z[I,0]<Z[I,1] * clusters are numbered from 0 to 2*NPoints-2, with indexes from 0 to NPoints-1 corresponding to points of the original dataset, and indexes from NPoints to 2*NPoints-2 correspond to clusters generated by subsequent merges (I-th row of Z creates cluster with index NPoints+I). IMPORTANT: indexes in Z[] are indexes in the ORIGINAL, unsorted dataset. In addition to Z algorithm outputs permutation which rearranges points in such way that subsequent merges are performed on adjacent points (such order is needed if you want to build dendrogram). However, indexes in Z are related to original, unrearranged sequence of points. P array[NPoints], permutation which reorders points for dendrogram construction. P[i] contains index of the position where we should move I-th point of the original dataset in order to apply merges PZ/PM. PZ same as Z, but for permutation of points given by P. The only thing which changed are indexes of the original points; indexes of clusters remained same. MergeDist array[NPoints-1], contains distances between clusters being merged (MergeDist[i] correspond to merge stored in Z[i,...]). PM array[NPoints-1,6], another representation of merges, which is suited for dendrogram construction. It deals with rearranged points (permutation P is applied) and represents merges in a form which different from one used by Z. For each I from 0 to NPoints-2, I-th row of PM represents merge performed on two clusters C0 and C1. Here: * C0 contains points with indexes PM[I,0]...PM[I,1] * C1 contains points with indexes PM[I,2]...PM[I,3] * indexes stored in PM are given for dataset sorted according to permutation P * PM[I,1]=PM[I,2]-1 (only adjacent clusters are merged) * PM[I,0]<=PM[I,1], PM[I,2]<=PM[I,3], i.e. both clusters contain at least one point * heights of "subdendrograms" corresponding to C0/C1 are stored in PM[I,4] and PM[I,5]. Subdendrograms corresponding to single-point clusters have height=0. Dendrogram of the merge result has height H=max(H0,H1)+1. NOTE: there is one-to-one correspondence between merges described by Z and PM. I-th row of Z describes same merge of clusters as I-th row of PM, with "left" cluster from Z corresponding to the "left" one from PM. -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public class ahcreport { public int npoints; public int[] p; public int[,] z; public int[,] pz; public int[,] pm; public double[] mergedist; }
clusterizerstate
class/************************************************************************* This structure is a clusterization engine. You should not try to access its fields directly. Use ALGLIB functions in order to work with this object. -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public class clusterizerstate { }
kmeansreport
class/************************************************************************* This structure is used to store results of the k-means++ clustering algorithm. Following information is always returned: * NPoints contains number of points in the original dataset * TerminationType contains completion code, negative on failure, positive on success * K contains number of clusters For positive TerminationType we return: * NFeatures contains number of variables in the original dataset * C, which contains centers found by algorithm * CIdx, which maps points of the original dataset to clusters FORMAL DESCRIPTION OF FIELDS: NPoints number of points, >=0 NFeatures number of variables, >=1 TerminationType completion code: * -5 if distance type is anything different from Euclidean metric * -3 for degenerate dataset: a) less than K distinct points, b) K=0 for non-empty dataset. * +1 for successful completion K number of clusters C array[K,NFeatures], rows of the array store centers CIdx array[NPoints], which contains cluster indexes -- ALGLIB -- Copyright 27.11.2012 by Bochkanov Sergey *************************************************************************/public class kmeansreport { public int npoints; public int nfeatures; public int terminationtype; public int k; public double[,] c; public int[] cidx; }
clusterizercreate
function/************************************************************************* This function initializes clusterizer object. Newly initialized object is empty, i.e. it does not contain dataset. You should use it as follows: 1. creation 2. dataset is added with ClusterizerSetPoints() 3. additional parameters are set 3. clusterization is performed with one of the clustering functions -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizercreate(out clusterizerstate s)
clusterizergetdistances
function/************************************************************************* This function returns distance matrix for dataset FOR USERS OF SMP EDITION: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Multicore version is pretty efficient on large ! problems which need more than 1.000.000 operations to be solved, ! gives moderate speed-up in mid-range (from 100.000 to 1.000.000 CPU ! cycles), but gives no speed-up for small problems (less than 100.000 ! operations). INPUT PARAMETERS: XY - array[NPoints,NFeatures], dataset NPoints - number of points, >=0 NFeatures- number of features, >=1 DistType- distance function: * 0 Chebyshev distance (L-inf norm) * 1 city block distance (L1 norm) * 2 Euclidean distance (L2 norm) * 10 Pearson correlation: dist(a,b) = 1-corr(a,b) * 11 Absolute Pearson correlation: dist(a,b) = 1-|corr(a,b)| * 12 Uncentered Pearson correlation (cosine of the angle): dist(a,b) = a'*b/(|a|*|b|) * 13 Absolute uncentered Pearson correlation dist(a,b) = |a'*b|/(|a|*|b|) * 20 Spearman rank correlation: dist(a,b) = 1-rankcorr(a,b) * 21 Absolute Spearman rank correlation dist(a,b) = 1-|rankcorr(a,b)| OUTPUT PARAMETERS: D - array[NPoints,NPoints], distance matrix (full matrix is returned, with lower and upper triangles) NOTES: different distance functions have different performance penalty: * Euclidean or Pearson correlation distances are the fastest ones * Spearman correlation distance function is a bit slower * city block and Chebyshev distances are order of magnitude slower The reason behing difference in performance is that correlation-based distance functions are computed using optimized linear algebra kernels, while Chebyshev and city block distance functions are computed using simple nested loops with two branches at each iteration. -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizergetdistances( double[,] xy, int npoints, int nfeatures, int disttype, out double[,] d) public static void alglib.smp_clusterizergetdistances( double[,] xy, int npoints, int nfeatures, int disttype, out double[,] d)
clusterizergetkclusters
function/************************************************************************* This function takes as input clusterization report Rep, desired clusters count K, and builds top K clusters from hierarchical clusterization tree. It returns assignment of points to clusters (array of cluster indexes). INPUT PARAMETERS: Rep - report from ClusterizerRunAHC() performed on XY K - desired number of clusters, 1<=K<=NPoints. K can be zero only when NPoints=0. OUTPUT PARAMETERS: CIdx - array[NPoints], I-th element contains cluster index (from 0 to K-1) for I-th point of the dataset. CZ - array[K]. This array allows to convert cluster indexes returned by this function to indexes used by Rep.Z. J-th cluster returned by this function corresponds to CZ[J]-th cluster stored in Rep.Z/PZ/PM. It is guaranteed that CZ[I]<CZ[I+1]. NOTE: K clusters built by this subroutine are assumed to have no hierarchy. Although they were obtained by manipulation with top K nodes of dendrogram (i.e. hierarchical decomposition of dataset), this function does not return information about hierarchy. Each of the clusters stand on its own. NOTE: Cluster indexes returned by this function does not correspond to indexes returned in Rep.Z/PZ/PM. Either you work with hierarchical representation of the dataset (dendrogram), or you work with "flat" representation returned by this function. Each of representations has its own clusters indexing system (former uses [0, 2*NPoints-2]), while latter uses [0..K-1]), although it is possible to perform conversion from one system to another by means of CZ array, returned by this function, which allows you to convert indexes stored in CIdx to the numeration system used by Rep.Z. NOTE: this subroutine is optimized for moderate values of K. Say, for K=5 it will perform many times faster than for K=100. Its worst-case performance is O(N*K), although in average case it perform better (up to O(N*log(K))). -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizergetkclusters( ahcreport rep, int k, out int[] cidx, out int[] cz)
clusterizerrunahc
function/************************************************************************* This function performs agglomerative hierarchical clustering FOR USERS OF SMP EDITION: ! This function can utilize multicore capabilities of your system. In ! order to do this you have to call version with "smp_" prefix, which ! indicates that multicore code will be used. ! !  We recommend you to carefully read ALGLIB Reference Manual, section ! called 'SMP support', before using parallel version of this function. ! ! You should remember that starting/stopping worker thread always have ! non-zero cost. Multicore version is pretty efficient on large ! problems which need more than 1.000.000 operations to be solved, ! gives moderate speed-up in mid-range (from 100.000 to 1.000.000 CPU ! cycles), but gives no speed-up for small problems (less than 100.000 ! operations). INPUT PARAMETERS: S - clusterizer state, initialized by ClusterizerCreate() OUTPUT PARAMETERS: Rep - clustering results; see description of AHCReport structure for more information. NOTE 1: hierarchical clustering algorithms require large amounts of memory. In particular, this implementation needs sizeof(double)*NPoints^2 bytes, which are used to store distance matrix. In case we work with user-supplied matrix, this amount is multiplied by 2 (we have to store original matrix and to work with its copy). For example, problem with 10000 points would require 800M of RAM, even when working in a 1-dimensional space. -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizerrunahc( clusterizerstate s, out ahcreport rep) public static void alglib.smp_clusterizerrunahc( clusterizerstate s, out ahcreport rep)
clusterizerrunkmeans
function/************************************************************************* This function performs clustering by k-means++ algorithm. You may change algorithm properties like number of restarts or iterations limit by calling ClusterizerSetKMeansLimits() functions. INPUT PARAMETERS: S - clusterizer state, initialized by ClusterizerCreate() K - number of clusters, K>=0. K can be zero only when algorithm is called for empty dataset, in this case completion code is set to success (+1). If K=0 and dataset size is non-zero, we can not meaningfully assign points to some center (there are no centers because K=0) and return -3 as completion code (failure). OUTPUT PARAMETERS: Rep - clustering results; see description of KMeansReport structure for more information. NOTE 1: k-means clustering can be performed only for datasets with Euclidean distance function. Algorithm will return negative completion code in Rep.TerminationType in case dataset was added to clusterizer with DistType other than Euclidean (or dataset was specified by distance matrix instead of explicitly given points). -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizerrunkmeans( clusterizerstate s, int k, out kmeansreport rep)
clusterizerseparatedbycorr
function/************************************************************************* This function accepts AHC report Rep, desired maximum intercluster correlation and returns top clusters from hierarchical clusterization tree which are separated by correlation R or LOWER. It returns assignment of points to clusters (array of cluster indexes). There is one more function with similar name - ClusterizerSeparatedByDist, which returns clusters with intercluster distance equal to R or HIGHER (note: higher for distance, lower for correlation). INPUT PARAMETERS: Rep - report from ClusterizerRunAHC() performed on XY R - desired maximum intercluster correlation, -1<=R<=+1 OUTPUT PARAMETERS: K - number of clusters, 1<=K<=NPoints CIdx - array[NPoints], I-th element contains cluster index (from 0 to K-1) for I-th point of the dataset. CZ - array[K]. This array allows to convert cluster indexes returned by this function to indexes used by Rep.Z. J-th cluster returned by this function corresponds to CZ[J]-th cluster stored in Rep.Z/PZ/PM. It is guaranteed that CZ[I]<CZ[I+1]. NOTE: K clusters built by this subroutine are assumed to have no hierarchy. Although they were obtained by manipulation with top K nodes of dendrogram (i.e. hierarchical decomposition of dataset), this function does not return information about hierarchy. Each of the clusters stand on its own. NOTE: Cluster indexes returned by this function does not correspond to indexes returned in Rep.Z/PZ/PM. Either you work with hierarchical representation of the dataset (dendrogram), or you work with "flat" representation returned by this function. Each of representations has its own clusters indexing system (former uses [0, 2*NPoints-2]), while latter uses [0..K-1]), although it is possible to perform conversion from one system to another by means of CZ array, returned by this function, which allows you to convert indexes stored in CIdx to the numeration system used by Rep.Z. NOTE: this subroutine is optimized for moderate values of K. Say, for K=5 it will perform many times faster than for K=100. Its worst-case performance is O(N*K), although in average case it perform better (up to O(N*log(K))). -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizerseparatedbycorr( ahcreport rep, double r, out int k, out int[] cidx, out int[] cz)
clusterizerseparatedbydist
function/************************************************************************* This function accepts AHC report Rep, desired minimum intercluster distance and returns top clusters from hierarchical clusterization tree which are separated by distance R or HIGHER. It returns assignment of points to clusters (array of cluster indexes). There is one more function with similar name - ClusterizerSeparatedByCorr, which returns clusters with intercluster correlation equal to R or LOWER (note: higher for distance, lower for correlation). INPUT PARAMETERS: Rep - report from ClusterizerRunAHC() performed on XY R - desired minimum intercluster distance, R>=0 OUTPUT PARAMETERS: K - number of clusters, 1<=K<=NPoints CIdx - array[NPoints], I-th element contains cluster index (from 0 to K-1) for I-th point of the dataset. CZ - array[K]. This array allows to convert cluster indexes returned by this function to indexes used by Rep.Z. J-th cluster returned by this function corresponds to CZ[J]-th cluster stored in Rep.Z/PZ/PM. It is guaranteed that CZ[I]<CZ[I+1]. NOTE: K clusters built by this subroutine are assumed to have no hierarchy. Although they were obtained by manipulation with top K nodes of dendrogram (i.e. hierarchical decomposition of dataset), this function does not return information about hierarchy. Each of the clusters stand on its own. NOTE: Cluster indexes returned by this function does not correspond to indexes returned in Rep.Z/PZ/PM. Either you work with hierarchical representation of the dataset (dendrogram), or you work with "flat" representation returned by this function. Each of representations has its own clusters indexing system (former uses [0, 2*NPoints-2]), while latter uses [0..K-1]), although it is possible to perform conversion from one system to another by means of CZ array, returned by this function, which allows you to convert indexes stored in CIdx to the numeration system used by Rep.Z. NOTE: this subroutine is optimized for moderate values of K. Say, for K=5 it will perform many times faster than for K=100. Its worst-case performance is O(N*K), although in average case it perform better (up to O(N*log(K))). -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizerseparatedbydist( ahcreport rep, double r, out int k, out int[] cidx, out int[] cz)
clusterizersetahcalgo
function/************************************************************************* This function sets agglomerative hierarchical clustering algorithm INPUT PARAMETERS: S - clusterizer state, initialized by ClusterizerCreate() Algo - algorithm type: * 0 complete linkage (default algorithm) * 1 single linkage * 2 unweighted average linkage * 3 weighted average linkage -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizersetahcalgo( clusterizerstate s, int algo)
clusterizersetdistances
function/************************************************************************* This function adds dataset given by distance matrix to the clusterizer structure. It is important that dataset is not given explicitly - only distance matrix is given. This function overrides all previous calls of ClusterizerSetPoints() or ClusterizerSetDistances(). INPUT PARAMETERS: S - clusterizer state, initialized by ClusterizerCreate() D - array[NPoints,NPoints], distance matrix given by its upper or lower triangle (main diagonal is ignored because its entries are expected to be zero). NPoints - number of points IsUpper - whether upper or lower triangle of D is given. NOTE 1: different clustering algorithms have different limitations: * agglomerative hierarchical clustering algorithms may be used with any kind of distance metric, including one which is given by distance matrix * k-means++ clustering algorithm may be used only with Euclidean distance function and explicitly given points - it can not be used with dataset given by distance matrix Thus, if you call this function, you will be unable to use k-means clustering algorithm to process your problem. -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizersetdistances( clusterizerstate s, double[,] d, bool isupper) public static void alglib.clusterizersetdistances( clusterizerstate s, double[,] d, int npoints, bool isupper)
Examples: [1]
clusterizersetkmeanslimits
function/************************************************************************* This function sets k-means++ properties : number of restarts and maximum number of iterations per one run. INPUT PARAMETERS: S - clusterizer state, initialized by ClusterizerCreate() Restarts- restarts count, >=1. k-means++ algorithm performs several restarts and chooses best set of centers (one with minimum squared distance). MaxIts - maximum number of k-means iterations performed during one run. >=0, zero value means that algorithm performs unlimited number of iterations. -- ALGLIB -- Copyright 10.07.2012 by Bochkanov Sergey *************************************************************************/public static void alglib.clusterizersetkmeanslimits( clusterizerstate s, int restarts, int maxits)
clusterizersetpoints