Module scemtools
source code
Implements the "Shuffled Complex Evolution Metropolis"
Algoritm of Vrugt et al. [1]
Reference:
[1] Jasper A. Vrugt, Hoshin V. Gupta, Willem Bouten, and Soroosh
Sorooshian A Shuffled Complex Evolution Metropolis algorithm for
optimization and uncertainty assessment of hydrologic model parameters,
WATER RESOURCES RESEARCH, VOL. 39, NO. 8, 1201, doi:10.1029/2002WR001642,
2003 http://www.agu.org/pubs/crossref/2003/2002WR001642.shtml
[2] Vrugt JA, Nuallain , Robinson BA, Bouten W, Dekker SC, Sloot PM
Application of parallel computing to stochastic parameter estimation in
environmental models, Computers & Geosciences, Vol. 32, No. 8.
(October 2006), pp. 11391155.
http://www.science.uva.nl/research/scs/papers/archive/Vrugt2006b.pdf

multinormal_pdf(mean,
var)
var must be symmetric positive definite 
source code



sequential_deal(inarray,
n)
inarray: should be a set of N objects (the objects can be vectors
themselves, but inarray should be indexable like a list. 
source code





sort_ab_with_b(a,
b,
ord=1)
default is descending... 
source code









sort_complex2(c,
a)
c and a are partially sorted (either the last one is bad, or the
first one) 
source code



update_complex(Ck,
ak,
c,
a,
pos)
ak is sorted (descending) 
source code



scem(Ck,
ak,
Sk,
Sak,
target,
cn)
This is the SCEM algorithm starting from line [35] of the reference [1]. 
source code





TYPE_DOUBLE = numpy.lib.getlimits.finfo(numpy.lib.getlimits.nu...


TINY = TYPE_DOUBLE.tiny


inarray: should be a set of N objects (the objects can be vectors
themselves, but inarray should be indexable like a list. It is
coerced into a numpy array because the last operations requires that
it is also indexable by a 'list.'

it should have a length divisble by n, otherwise the reshape will
fail (this is a feature !)

sequential_deal(range(20), 5) wil return a 5 element list, each
element being a 4list of index. (see below)
>>> for l in sequential_deal(range(20),5):
... print l
...
[ 0 5 10 15]
[ 1 6 11 16]
[ 2 7 12 17]
[ 3 8 13 18]
[ 4 9 14 19]

This is the SCEM algorithm starting from line [35] of the reference [1].
 [inout] Ck is the kth 'complex' with m points. This should be an m by n array
n being the dimensionality of the density function. i.e., the data
are arranged in rows.
Ck is assumed to be sorted according to the target density.
 [inout] ak, the density of the points of Ck.
 [inout] Sk, the entire chain. (should be a list)
 [inout] Sak, the cost of the entire chain (should be a list)
Sak would be more convenient to use if it is a numpy array, but we need
to append to it frequently.
 [in] target: target density function
 [in] cn: jumprate. (see Paragraph 37 of [1.]
 The invariants: ak is always aligned with Ck, and are the cost of Ck
 Similarly, Sak is always aligned with Sk in the same way.
 On return... sort order in Ck/ak is destroyed. but see sort_complex2

Mixing and dealing the complexes. The types of Cs and As are very
important....

TYPE_DOUBLE
 Value:
numpy.lib.getlimits.finfo(numpy.lib.getlimits.numeric.double)

