This webpage is for programmers who need examples of use of the functions of the module. The examples are designed to illustrate the syntax. They do not correspond to any meaningful model. For examples of models, visit biogeme.epfl.ch.
import datetime
print(datetime.datetime.now())
2021-08-04 16:38:09.707948
import biogeme.version as ver
print(ver.getText())
biogeme 3.2.8 [2021-08-04] Version entirely written in Python Home page: http://biogeme.epfl.ch Submit questions to https://groups.google.com/d/forum/biogeme Michel Bierlaire, Transport and Mobility Laboratory, Ecole Polytechnique Fédérale de Lausanne (EPFL)
import numpy as np
We set the seed so that the outcome of random operations is always the same.
np.random.seed(90267)
import biogeme.draws as dr
draws = dr.getUniform(sampleSize=3,
numberOfDraws=10,
symmetric=False)
draws
array([[0.4069955 , 0.95209084, 0.39330002, 0.42242449, 0.39231397, 0.65602235, 0.07182417, 0.28215248, 0.29462105, 0.29437835], [0.71971885, 0.97647083, 0.86001524, 0.63780397, 0.09774763, 0.21585889, 0.90199693, 0.25945665, 0.24673199, 0.60813096], [0.9140627 , 0.25344881, 0.29963774, 0.02324269, 0.00851069, 0.3091653 , 0.17208375, 0.84805301, 0.77527991, 0.23414075]])
draws = dr.getUniform(sampleSize=3,
numberOfDraws=10,
symmetric=True)
draws
array([[ 0.54159507, -0.61343854, 0.41089924, 0.65656456, -0.37767775, 0.9788048 , -0.19924935, -0.49970433, -0.94157778, 0.43328274], [-0.39595535, 0.44247489, 0.26526338, 0.51113532, -0.50745861, 0.73837224, -0.86845256, -0.62098784, 0.55355381, -0.21390505], [-0.7972486 , 0.43932172, 0.12704909, 0.97976216, 0.60990425, -0.43047078, 0.3801898 , 0.53669883, 0.02269837, -0.76152749]])
The Modified Latin Hypercube Sampling (MLHS, Hess et al, 2006) provides U[0,1] draws from a perturbed grid, designed for Monte-Carlo integration.
latinHypercube = dr.getLatinHypercubeDraws(sampleSize=3,
numberOfDraws=10)
latinHypercube
array([[0.32521073, 0.63000128, 0.17746161, 0.53094562, 0.41705476, 0.82583115, 0.78232205, 0.08874858, 0.15503933, 0.53946678], [0.9864292 , 0.02833746, 0.58581152, 0.84659261, 0.66113494, 0.92572347, 0.39798632, 0.35041124, 0.9553105 , 0.70211636], [0.74494216, 0.88475329, 0.04655552, 0.13146474, 0.46312566, 0.28911263, 0.68569908, 0.25111488, 0.204371 , 0.48178854]])
The same method can be used to generate draws from U[-1,1]
latinHypercube = dr.getLatinHypercubeDraws(sampleSize=5,
numberOfDraws=10,
symmetric=True)
latinHypercube
array([[-0.48782677, 0.94658307, -0.29089021, -0.99773248, 0.68546628, -0.92901557, 0.81537697, 0.06090223, -0.70370281, -0.66028972], [-0.11977072, 0.90385423, 0.5185657 , -0.17531685, -0.59367191, -0.4766856 , -0.24987929, 0.01404931, 0.59241613, 0.33354383], [-0.53646676, -0.60271978, 0.42173134, 0.25802578, -0.36910209, -0.21567312, 0.9601231 , -0.40689483, 0.65280404, 0.18610977], [ 0.14925908, -0.0534143 , -0.35650472, 0.63648376, 0.52128198, -0.86688846, -0.78450064, 0.31241634, -0.74722651, 0.11262522], [-0.91420031, 0.23854628, 0.73852358, -0.02715404, 0.79892765, -0.13414034, -0.80093917, 0.46081992, 0.36229167, 0.85247066]])
The user can provide her own series of U[0,1] draws.
myUnif = np.random.uniform(size=30)
myUnif
array([0.06175992, 0.20160719, 0.77983176, 0.6025268 , 0.696678 , 0.76427402, 0.56865543, 0.37119385, 0.23323665, 0.0236732 , 0.33989596, 0.06099459, 0.2213469 , 0.47601499, 0.86605521, 0.35092633, 0.64428447, 0.90464931, 0.85755345, 0.13233163, 0.4088656 , 0.41519067, 0.50233919, 0.49837905, 0.61238122, 0.67505865, 0.34980957, 0.03300237, 0.14517637, 0.23786104])
latinHypercube = dr.getLatinHypercubeDraws(sampleSize=3,
numberOfDraws=10,
symmetric=False,
uniformNumbers=myUnif)
latinHypercube
array([[0.7832793 , 0.09266106, 0.40737823, 0.27444122, 0.62858511, 0.90110008, 0.59682164, 0.30078911, 0.24570646, 0.21895518], [0.00205866, 0.75007797, 0.55480948, 0.15655593, 0.97459537, 0.19214247, 0.82041271, 0.04005357, 0.68029552, 0.93817255], [0.51169754, 0.4492005 , 0.63774439, 0.12008423, 0.85583529, 0.87832699, 0.36869982, 0.71383969, 0.3446632 , 0.49553517]])
The uniform draws can also be arranged in a two-dimension array
myUnif = dr.getUniform(sampleSize=3,
numberOfDraws=10)
myUnif
array([[0.42934713, 0.60387549, 0.8965652 , 0.75460304, 0.76976185, 0.45135228, 0.83200157, 0.22472631, 0.90449758, 0.09730533], [0.09163669, 0.28808628, 0.13395278, 0.36000124, 0.79323743, 0.56924873, 0.33890064, 0.98226408, 0.66005216, 0.07526522], [0.33633631, 0.37411379, 0.46400897, 0.22787268, 0.0081577 , 0.28662478, 0.93417098, 0.10990211, 0.94595662, 0.75624039]])
latinHypercube = dr.getLatinHypercubeDraws(sampleSize=3,
numberOfDraws=10,
uniformNumbers=myUnif)
latinHypercube
array([[0.18171174, 0.80027192, 0.9036634 , 0.8978057 , 0.49310791, 0.51897496, 0.77426242, 0.54463002, 0.5994088 , 0.12515343], [0.67787788, 0.99187468, 0.7488003 , 0.33638789, 0.44533337, 0.40446509, 0.71247046, 0.15899206, 0.01431157, 0.63584217], [0.09655217, 0.29681659, 0.05346252, 0.24082421, 0.30324351, 0.84288749, 0.22773339, 0.37626954, 0.96486522, 0.62200174]])
One Halton sequence
halton = dr.getHaltonDraws(sampleSize=2,
numberOfDraws=10,
base=3)
halton
array([[0.33333333, 0.66666667, 0.11111111, 0.44444444, 0.77777778, 0.22222222, 0.55555556, 0.88888889, 0.03703704, 0.37037037], [0.7037037 , 0.14814815, 0.48148148, 0.81481481, 0.25925926, 0.59259259, 0.92592593, 0.07407407, 0.40740741, 0.74074074]])
Several Halton sequences
halton = dr.getHaltonDraws(sampleSize=3,
numberOfDraws=10)
halton
array([[0.5 , 0.25 , 0.75 , 0.125 , 0.625 , 0.375 , 0.875 , 0.0625 , 0.5625 , 0.3125 ], [0.8125 , 0.1875 , 0.6875 , 0.4375 , 0.9375 , 0.03125, 0.53125, 0.28125, 0.78125, 0.15625], [0.65625, 0.40625, 0.90625, 0.09375, 0.59375, 0.34375, 0.84375, 0.21875, 0.71875, 0.46875]])
Shuffled Halton sequences
halton = dr.getHaltonDraws(sampleSize=3,
numberOfDraws=10,
shuffled=True)
halton
array([[0.25 , 0.125 , 0.28125, 0.46875, 0.53125, 0.59375, 0.09375, 0.3125 , 0.8125 , 0.21875], [0.71875, 0.78125, 0.5 , 0.375 , 0.9375 , 0.625 , 0.875 , 0.90625, 0.40625, 0.65625], [0.84375, 0.1875 , 0.4375 , 0.5625 , 0.34375, 0.03125, 0.6875 , 0.15625, 0.0625 , 0.75 ]])
The above sequences were generated using the default base: 2. It is possible to generate sequences using different prime numbers.
halton = dr.getHaltonDraws(sampleSize=1,
numberOfDraws=10,
base=3)
halton
array([[0.33333333, 0.66666667, 0.11111111, 0.44444444, 0.77777778, 0.22222222, 0.55555556, 0.88888889, 0.03703704, 0.37037037]])
It is also possible to skip the first items of the sequence. This is desirable in the context of Monte-Carlo integration.
halton = dr.getHaltonDraws(sampleSize=1,
numberOfDraws=10,
base=3,
skip=10)
halton
array([[0.7037037 , 0.14814815, 0.48148148, 0.81481481, 0.25925926, 0.59259259, 0.92592593, 0.07407407, 0.40740741, 0.74074074]])
Antithetic draws can be generated from any function generating uniform draws
draws = dr.getAntithetic(dr.getUniform,
sampleSize=3,
numberOfDraws=10)
draws
array([[0.42212598, 0.47179374, 0.148715 , 0.73860092, 0.04105165, 0.57787402, 0.52820626, 0.851285 , 0.26139908, 0.95894835], [0.51232324, 0.76571874, 0.50131744, 0.1464681 , 0.55853883, 0.48767676, 0.23428126, 0.49868256, 0.8535319 , 0.44146117], [0.41770545, 0.77892061, 0.75202037, 0.04320253, 0.25552776, 0.58229455, 0.22107939, 0.24797963, 0.95679747, 0.74447224]])
Antithetic MLHS
draws = dr.getAntithetic(dr.getLatinHypercubeDraws,
sampleSize=3,
numberOfDraws=10)
draws
array([[0.0780733 , 0.31379531, 0.15866299, 0.0151855 , 0.42181169, 0.9219267 , 0.68620469, 0.84133701, 0.9848145 , 0.57818831], [0.47828697, 0.76997285, 0.24430161, 0.61220294, 0.90674633, 0.52171303, 0.23002715, 0.75569839, 0.38779706, 0.09325367], [0.94935667, 0.56168759, 0.72254999, 0.35044001, 0.80752534, 0.05064333, 0.43831241, 0.27745001, 0.64955999, 0.19247466]])
Antithetic Halton.
draws = dr.getAntithetic(dr.getHaltonDraws,
sampleSize=1,
numberOfDraws=10)
draws
array([[0.5 , 0.25 , 0.75 , 0.125, 0.625, 0.5 , 0.75 , 0.25 , 0.875, 0.375]])
As antithetic Halton draws may be correlated, it is a good idea to skip the first draws
def halton(sampleSize,numberOfDraws):
return dr.getHaltonDraws(numberOfDraws,
sampleSize,skip=100)
draws = dr.getAntithetic(halton,
sampleSize=3,
numberOfDraws=10)
draws
array([[0.6484375, 0.3984375, 0.8984375, 0.3515625, 0.6015625, 0.1015625], [0.0859375, 0.5859375, 0.3359375, 0.9140625, 0.4140625, 0.6640625], [0.8359375, 0.2109375, 0.7109375, 0.1640625, 0.7890625, 0.2890625], [0.4609375, 0.9609375, 0.0546875, 0.5390625, 0.0390625, 0.9453125], [0.5546875, 0.3046875, 0.8046875, 0.4453125, 0.6953125, 0.1953125]])
Generate pseudo-random numbers from a normal distribution N(0,1) using the Algorithm AS241 Appl. Statist. (1988) Vol. 37, No. 3 by Wichura
draws = dr.getNormalWichuraDraws(sampleSize=3,
numberOfDraws=10)
draws
array([[ 0.46838353, 0.30944461, -0.72705093, 1.47731844, 2.32498619, -0.38523969, 2.41753668, -0.89851703, 1.04351099, -0.52799835], [ 1.92824184, 0.39265187, -1.69805707, 0.72243106, 0.23153714, 1.29178853, -0.39065241, 0.75741968, -0.1117659 , -0.25373762], [-1.62848672, 0.80968139, -0.47731943, -0.6055917 , -1.04245527, 2.17223027, -1.20160168, 1.83144122, -1.49289716, 0.7304994 ]])
The antithetic version actually generates half of the draws and complete them with their antithetic version
draws = dr.getNormalWichuraDraws(sampleSize=3,
numberOfDraws=10,
antithetic=True)
draws
array([[-0.39745017, -0.39736937, 1.09672047, -0.51651723, 0.25290974, 0.39745017, 0.39736937, -1.09672047, 0.51651723, -0.25290974], [-0.3183153 , 0.48156301, 0.99268525, 0.36495778, 0.60162053, 0.3183153 , -0.48156301, -0.99268525, -0.36495778, -0.60162053], [-0.08619037, 0.78059876, -0.17844824, 0.16459028, -0.03529156, 0.08619037, -0.78059876, 0.17844824, -0.16459028, 0.03529156]])
The user can provide her own series of U[0,1] draws. In this example, we use the MLHS procedure to generate these draws. Note that, if the antithetic version is used, only half of the requested draws must be provided.
myUnif = dr.getLatinHypercubeDraws(sampleSize=3,
numberOfDraws=5)
myUnif
array([[0.73078091, 0.97414181, 0.54479742, 0.0906173 , 0.65406372], [0.52895577, 0.22199862, 0.00909876, 0.42051145, 0.87765406], [0.15495423, 0.84811817, 0.76065654, 0.33554539, 0.32617466]])
draws = dr.getNormalWichuraDraws(sampleSize=3,
numberOfDraws=10,
uniformNumbers=myUnif,
antithetic=True)
draws
array([[ 0.61517648, 1.94548717, 0.11252751, -1.33696337, 0.39631515, -0.61517648, -1.94548717, -0.11252751, 1.33696337, -0.39631515], [ 0.0726452 , -0.76546072, -2.36112653, -0.20058523, 1.16333929, -0.0726452 , 0.76546072, 2.36112653, 0.20058523, -1.16333929], [-1.01541415, 1.02839586, 0.70841608, -0.42465146, -0.45050087, 1.01541415, -1.02839586, -0.70841608, 0.42465146, 0.45050087]])
The same with Halton draws
myUnif = dr.getHaltonDraws(sampleSize=2,
numberOfDraws=5,
base=3,
skip=10)
myUnif
array([[0.7037037 , 0.14814815, 0.48148148, 0.81481481, 0.25925926], [0.59259259, 0.92592593, 0.07407407, 0.40740741, 0.74074074]])
draws = dr.getNormalWichuraDraws(numberOfDraws=10,
sampleSize=2,
uniformNumbers=myUnif,
antithetic=True)
draws
array([[ 0.53508282, -1.04440879, -0.04643574, 0.89577982, -0.64563075, -0.53508282, 1.04440879, 0.04643574, -0.89577982, 0.64563075], [ 0.2342192 , 1.44610359, -1.44610359, -0.23421919, 0.64563075, -0.2342192 , -1.44610359, 1.44610359, 0.23421919, -0.64563075]])