def test_sens(self):
"""Test the sign and magnitude of the sensitivity
"""
initial_data = {'simple': self.sim1({'simp': self.model})}
sens_matrix = self.bayes._get_sens(initial_data=initial_data)
indep = np.arange(10)
resp_mat = np.array([(1.02 * 2 * indep)**2 + 1 * indep -
(2 * indep)**2 - indep,
(2 * indep)**2 + 1.02 * indep -
(2 * indep)**2 - 1 * indep])
inp_mat = np.array([[0.02 * 2, 0], [0, 0.02]])
true_sens = np.linalg.lstsq(inp_mat, resp_mat)[0].T
npt.assert_array_almost_equal(sens_matrix['simple'], true_sens,
decimal=8,
err_msg='Sens matrix not as predicted')
resp_mat2 = resp_mat
resp_mat2[0, :] /= 0.04
resp_mat2[1, :] /= 0.02
npt.assert_array_almost_equal(sens_matrix['simple'], resp_mat2.T,
decimal=8,
err_msg='Sens matrix not as predicted')
python类assert_array_almost_equal()的实例源码
def test_sens_passed_inp_data(self):
"""Test the sensitivity using provided initial data
"""
data = self.bayes.get_data()
sens_matrix = self.bayes._get_sens(initial_data=data)
indep = np.linspace(1, 11, 5) - data['simple'][2]['tau']
resp_mat = np.array([(1.02 * 2 * indep)**2 + 1 * indep -
(2 * indep)**2 - indep,
(2 * indep)**2 + 1.02 * indep -
(2 * indep)**2 - 1 * indep])
inp_mat = np.array([[0.02 * 2, 0], [0, 0.02]])
true_sens = np.linalg.lstsq(inp_mat, resp_mat)[0].T
npt.assert_array_almost_equal(sens_matrix['simple'], true_sens,
decimal=8,
err_msg='Sens matrix not as predicted')
def test_model_pq(self):
"""Test the model PQ matrix generation
"""
new_model = self.bayes.update(models={
'simp': self.model.update_dof([4, 2])})
P, q = new_model._get_model_pq()
epsilon = np.array([2, 1])
sigma = inv(np.diag(np.ones(2)))
P_true = sigma
q_true = -np.dot(epsilon, sigma)
npt.assert_array_almost_equal(P, P_true, decimal=8)
npt.assert_array_almost_equal(q, q_true, decimal=8)
def test_sim_pq(self):
"""Test the simulation PQ matrix generation
"""
initial_data = self.bayes.get_data()
sens_matrix = self.bayes._get_sens(initial_data=initial_data)
P, q = self.bayes._get_sim_pq(initial_data, sens_matrix)
sigma = self.exp1.get_sigma()
P_true = np.dot(np.dot(sens_matrix['simple'].T,
inv(sigma)),
sens_matrix['simple'])
npt.assert_array_almost_equal(P, P_true, decimal=8)
epsilon = self.bayes.simulations['simple']['exp']\
.compare(initial_data['simple'])
q_true = -np.dot(np.dot(epsilon, inv(sigma)), sens_matrix['simple'])
npt.assert_array_almost_equal(q, q_true, decimal=8)
def test_compare(self):
"""Test the compare function
"""
xx = np.arange(10)
yy = (3 * xx)**2 + xx
models = {'simp': SimpleModel([2, 1])}
model_data = self.simSimp(models)
print(model_data)
print(yy)
res = self.simSimp.compare(model_data, (xx, (yy,), {'tau': 0}))
self.assertIsInstance(res, np.ndarray)
self.assertEqual(len(res), 10)
npt.assert_array_almost_equal(5 * xx**2, res)
def test_pll_sens(self):
"""Test sensitivity calculation
"""
models = {'simp': SimpleModel([2, 1])}
t0 = time.time()
sens = self.simSimp.get_sens_pll(models, 'simp')
t0 = time.time() - t0
print('Parallel sense took {:f}'.format(t0))
self.assertIsInstance(sens, np.ndarray)
self.assertEqual(sens.shape, (10, 2))
indep = np.arange(10)
resp_mat = np.array([(1.02 * 2 * indep)**2 + 1 * indep -
(2 * indep)**2 - indep,
(2 * indep)**2 + 1.02 * indep -
(2 * indep)**2 - 1 * indep])
inp_mat = np.array([[0.02 * 2, 0], [0, 0.02]])
true_sens = np.linalg.lstsq(inp_mat, resp_mat)[0].T
npt.assert_array_almost_equal(sens, true_sens, decimal=8)
def test_rotate_certain_angle_2():
d1 = np.array([1.0, 0.0])
d2 = np.array([0.0, 1.0])
dnew1, dnew2 = rotate.rotate_certain_angle(d1, d2, 90.0)
npt.assert_array_almost_equal(dnew1, [0.0, 1.0])
npt.assert_array_almost_equal(dnew2, [-1.0, 0.0])
dnew1, dnew2 = rotate.rotate_certain_angle(d1, d2, 180.0)
npt.assert_array_almost_equal(dnew1, [-1.0, 0.0])
npt.assert_array_almost_equal(dnew2, [0.0, -1.0])
dnew1, dnew2 = rotate.rotate_certain_angle(d1, d2, 270.0)
npt.assert_array_almost_equal(dnew1, [0.0, -1.0])
npt.assert_array_almost_equal(dnew2, [1.0, 0.0])
dnew1, dnew2 = rotate.rotate_certain_angle(d1, d2, 360.0)
npt.assert_array_almost_equal(dnew1, [1.0, 0.0])
npt.assert_array_almost_equal(dnew2, [0.0, 1.0])
def test_dump_adjsrc():
array = np.array([1., 2., 3., 4., 5.])
adj = AdjointSource(
"cc_traveltime_misfit", 2.0, 1.0, 17, 40, "BHZ",
adjoint_source=array, network="II", station="AAK",
location="", starttime=UTCDateTime(1990, 1, 1))
station_info = {"latitude": 1.0, "longitude": 2.0, "depth_in_m": 3.0,
"elevation_in_m": 4.0}
adj_array, adj_path, parameters = sa.dump_adjsrc(adj, station_info)
npt.assert_array_almost_equal(adj_array, array)
for key in station_info:
npt.assert_almost_equal(station_info[key], parameters[key])
assert adj_path == "II_AAK_BHZ"
npt.assert_almost_equal(parameters["misfit"], 2.0)
npt.assert_almost_equal(parameters["dt"], 1.0)
npt.assert_almost_equal(parameters["min_period"], 17.0)
npt.assert_almost_equal(parameters["max_period"], 40.0)
assert parameters["adjoint_source_type"] == "cc_traveltime_misfit"
assert parameters["station_id"] == "II.AAK"
assert parameters["component"], "BHZ"
assert UTCDateTime(parameters["starttime"]) == UTCDateTime(1990, 1, 1)
assert parameters["units"] == "m"
def test_cal_price_fft(self):
fft_pricer = FourierPricer(self.vanilla_option)
strike_arr = np.array([5, 10, 30, 36, 50, 60, 100])
put_call_arr = np.array(['call', 'call', 'put', 'call', 'call', 'put', 'call'])
exp = np.array(
[3.09958567e+01, 2.60163625e+01, 8.25753140e-05, 8.12953226e-01,
8.97449491e-11, 2.37785797e+01, 2.19293560e-85, ]
)
volatility = 0.20
N = 2 ** 15
d_u = 0.01
alpha = 1
fft_pricer.set_log_st_process(BlackScholes(volatility))
fft_pricer.set_pricing_engine(FFTEngine(N, d_u, alpha, spline_order=3))
res = fft_pricer.calc_price(strike_arr, put_call_arr, put_label='put')
npt.assert_array_almost_equal(res, exp, 6)
def test_cal_price_fractional_fft(self):
fft_pricer = FourierPricer(self.vanilla_option)
strike_arr = np.array([5, 10, 30, 36, 50, 60, 100])
put_call_arr = np.array(['call', 'call', 'put', 'call', 'call', 'put', 'call'])
exp = np.array(
[3.09958567e+01, 2.60163625e+01, 8.25753140e-05, 8.12953226e-01,
8.97449491e-11, 2.37785797e+01, 2.19293560e-85, ]
)
volatility = 0.20
N = 2 ** 14
d_u = 0.01
d_k = 0.01
alpha = 1
fft_pricer.set_log_st_process(BlackScholes(volatility))
fft_pricer.set_pricing_engine(FractionFFTEngine(N, d_u, d_k, alpha, spline_order=3))
res = fft_pricer.calc_price(strike_arr, put_call_arr, put_label='put')
npt.assert_array_almost_equal(res, exp, 6)
def test_cal_price_cosine_2(self):
cosine_pricer = FourierPricer(self.vanilla_option)
strike_arr = np.array([5, 10, 30, 36, 50, 60])
put_call_arr = np.array(['call', 'call', 'put', 'call', 'call', 'put'])
exp = np.array([3.09958568e+01, 2.60164423e+01, 9.56077953e-02, 8.81357807e-01,
1.41769466e-10, 2.37785797e+01]
)
volatility = 0.20
N = 2000
cosine_pricer.set_log_st_process(MertonJump(sigma=volatility,
jump_rate=0.090913148257155449,
norm_m=-0.91157356544103341,
norm_sig=7.3383200797618833e-05))
cosine_pricer.set_pricing_engine(CosineEngine(N, L=30))
res = cosine_pricer.calc_price(strike_arr, put_call_arr, put_label='put')
npt.assert_array_almost_equal(res, exp, 6)
def test_cal_price_cosine_3(self):
cosine_pricer = FourierPricer(self.vanilla_option)
strike_arr = np.array([5, 10, 30, 36, 40, 60])
put_call_arr = np.array(['call', 'call', 'put', 'call', 'call', 'put'])
exp = np.array([3.09958567e+01, 2.60163625e+01, 1.71886506e-04, 8.75203272e-01,
3.55292239e-02, 2.37785797e+01]
)
volatility = 0.20
N = 2000
cosine_pricer.set_log_st_process(KouJump(sigma=volatility,
jump_rate=23.339325557373201,
exp_pos=59.378410421004197,
exp_neg=-59.447921334340137,
prob_pos=-200.08018971817182))
cosine_pricer.set_pricing_engine(CosineEngine(N, L=30))
res = cosine_pricer.calc_price(strike_arr, put_call_arr, put_label='put')
npt.assert_array_almost_equal(res, exp, 6)
def test_cal_price_cosine_5(self):
cosine_pricer = FourierPricer(self.vanilla_option)
strike_arr = np.array([5, 10, 30, 36, 40, 60])
put_call_arr = np.array(['call', 'call', 'put', 'call', 'call', 'put'])
exp = np.array([22.48662613, 17.50712233, -8.50824394, -7.13267131, -7.22087064,
16.09965887]
)
volatility = 0.20
N = 2000
cosine_pricer.set_log_st_process(Poisson(
jump_rate=0.339325557373201,
))
cosine_pricer.set_pricing_engine(CosineEngine(N, L=30))
res = cosine_pricer.calc_price(strike_arr, put_call_arr, put_label='put')
npt.assert_array_almost_equal(res, exp, 6)
def test_cal_price_cosine_6(self):
cosine_pricer = FourierPricer(self.vanilla_option)
strike_arr = np.array([5, 10, 30, 36, 40, 60])
put_call_arr = np.array(['call', 'call', 'put', 'call', 'call', 'put'])
exp = np.array([30.9958673, 26.01656399, 0.15928293, 1.72971868, 0.56756891, 23.82357896])
volatility = 0.20
N = 2000
cosine_pricer.set_log_st_process(CGMY(c=0.1,
g=2,
m=2,
y=1.5))
cosine_pricer.set_pricing_engine(CosineEngine(N, L=30))
res = cosine_pricer.calc_price(strike_arr, put_call_arr, put_label='put')
npt.assert_array_almost_equal(res, exp, 6)
def compReadWrite(self, testMage, casttype=None, compressor=None, clevel = 1 ):
# This is the main functions which reads and writes from disk.
mrcName = os.path.join( tmpDir, "testMage.mrc" )
pixelsize = np.array( [1.2, 2.6, 3.4] )
mrcz.writeMRC( testMage, mrcName, dtype=casttype,
pixelsize=pixelsize, pixelunits=u"\AA",
voltage=300.0, C3=2.7, gain=1.05,
compressor=compressor, clevel=clevel )
rereadMage, rereadHeader = mrcz.readMRC( mrcName, pixelunits=u"\AA")
try: os.remove( mrcName )
except IOError: log.info( "Warning: file {} left on disk".format(mrcName) )
npt.assert_array_almost_equal( testMage, rereadMage )
npt.assert_array_equal( rereadHeader['voltage'], 300.0 )
npt.assert_array_almost_equal( rereadHeader['pixelsize'], pixelsize )
npt.assert_array_equal( rereadHeader['pixelunits'], u"\AA" )
npt.assert_array_equal( rereadHeader['C3'], 2.7 )
npt.assert_array_equal( rereadHeader['gain'], 1.05 )
def test_sigmoid_activation():
x2 = ad.Variable(name='x2')
y = au.nn.sigmoid(x2)
x2_val = np.array([-100, 0, 100])
grad_x2, = ad.gradients(y, [x2])
executor = ad.Executor([y, grad_x2])
y_val, grad_x2_val = executor.run(feed_shapes={x2: x2_val})
npt.assert_array_almost_equal(np.array([0.000, 0.500, 1.0]), y_val)
npt.assert_array_almost_equal(np.array([0, 0.25, 0]), grad_x2_val)
# testing with extreme values for numerical stability.
x2_val = np.array([-9.9e10, 9.9e10]).astype(np.float32)
y_val, grad_x2_val = executor.run(feed_shapes={x2: x2_val})
npt.assert_array_almost_equal(np.array([0.0, 1.0]), y_val)
npt.assert_array_almost_equal(np.array([0.0, 0.0]), grad_x2_val)
def test_max_pooling():
x2 = ad.Variable(name='x2')
y = au.nn.maxPool(x2, filter=(2, 2), strides=(2, 2))
grad_x2, = ad.gradients(y, [x2])
executor = ad.Executor([y, grad_x2])
x2_val = np.random.randn(1, 1, 4, 4)
y_val, grad_x2_val = executor.run(feed_shapes={x2: x2_val})
numerical_grad_x2 = ad.eval_numerical_grad(y,
feed_dict={x2: x2_val},
wrt=x2,
h=1e-5)
assert isinstance(y, ad.Node)
# TODO: (upul) looks like a bug in my eval_numerical_grad implementation
# Hence I'm using one decimal points
npt.assert_array_almost_equal(grad_x2_val, numerical_grad_x2, decimal=1)
def assert_array_collections_equal(correct, test, decimal=7):
"""Assert that two collections of numpy arrays have the same values.
Collections can be either a Sequence or a Mapping.
"""
if type(correct) != type(test):
raise ValueError('correct ({}) and test ({}) must have the same type.'.format(type(correct), type(test)))
assert_equal = lambda c, t: assert_array_almost_equal(c, t, decimal=decimal)
if isinstance(correct, Sequence):
assert len(correct) == len(test)
for c, t in izip(correct, test):
assert_equal(c, t)
elif isinstance(correct, Mapping):
# same keys
assert set(test.keys()) == set(correct.keys())
# same values
for key in test:
assert_equal(correct[key], test[key])
else:
raise TypeError('Inputs must be of type Mapping or Sequence, not {}.'.format(type(correct)))
def assert_array_collections_equal(correct, test, decimal=7):
"""Assert that two collections of numpy arrays have the same values.
Collections can be either a Sequence or a Mapping.
"""
if type(correct) != type(test):
raise ValueError('correct ({}) and test ({}) must have the same type.'.format(type(correct), type(test)))
assert_equal = lambda c, t: assert_array_almost_equal(c, t, decimal=decimal)
if isinstance(correct, Sequence):
assert len(correct) == len(test)
for c, t in izip(correct, test):
assert_equal(c, t)
elif isinstance(correct, Mapping):
# same keys
assert set(test.keys()) == set(correct.keys())
# same values
for key in test:
assert_equal(correct[key], test[key])
else:
raise TypeError('Inputs must be of type Mapping or Sequence, not {}.'.format(type(correct)))
def test_beta(self):
"""check that numerical output of Python code equals Mathematica code"""
my_beta = beta.beta(C, HIGHSCALE)
self.assertEqual(list(my_beta.keys()), beta.C_keys)
for k in beta.C_keys:
if isinstance(my_beta[k], float) or isinstance(my_beta[k], complex):
self.assertEqual(beta.C_keys_shape[k], 1)
else:
self.assertEqual(my_beta[k].shape, beta.C_keys_shape[k], msg=k)
for i, n in enumerate(C0):
self.assertAlmostEqual(betas_re[0][i]/my_beta[n].real, 1, places=4)
for i, n in enumerate(C2):
npt.assert_array_almost_equal(betas_re[1][i]/my_beta[n].real, np.ones((3,3)), decimal=2)
for i, n in enumerate(C4):
npt.assert_array_almost_equal(betas_re[2][i]/my_beta[n].real, np.ones((3,3,3,3)), decimal=2)
for i, n in enumerate(C2):
npt.assert_array_almost_equal(betas_im[1][i]/my_beta[n].imag, np.ones((3,3)), decimal=2)
for i, n in enumerate(C4):
npt.assert_array_almost_equal(betas_im[2][i]/my_beta[n].imag, np.ones((3,3,3,3)), decimal=2)
def test_predict_proba_rocauc(self):
"""
Test ROCAUC with classifiers that utilize predict_proba
"""
# Load the model and assert there is no decision_function method.
model = MultinomialNB()
with self.assertRaises(AttributeError):
model.decision_function
# Fit model and visualizer
visualizer = ROCAUC(model)
visualizer.fit(X, yb)
expected = np.asarray([
[0.493788, 0.506212],
[0.493375, 0.506625],
[0.493572, 0.506428],
[0.511063, 0.488936],
[0.511887, 0.488112],
[0.510841, 0.489158],
])
# Get the predict_proba scores and evaluate
y_scores = visualizer._get_y_scores(X)
npt.assert_array_almost_equal(y_scores, expected)
def test_download_prediction_csv_regr(driver, project, dataset, featureset, model, prediction):
driver.get('/')
_click_download(project.id, driver)
assert os.path.exists('/tmp/cesium_prediction_results.csv')
try:
results = np.genfromtxt('/tmp/cesium_prediction_results.csv',
dtype='str', delimiter=',')
npt.assert_equal(results[0],
['ts_name', 'label', 'prediction'])
npt.assert_array_almost_equal(
[[float(e) for e in row] for row in results[1:]],
[[0, 2.2, 2.2],
[1, 3.4, 3.4],
[2, 4.4, 4.4],
[3, 2.2, 2.2],
[4, 3.1, 3.1]])
finally:
os.remove('/tmp/cesium_prediction_results.csv')
def compare_cwt():
"""Compare the output of Scipy's cwt (using direct convolution)
and my cwt (using fft convolution).
"""
cwt = scipy.signal.cwt
fft_cwt = wavelets.cwt
data = np.random.random(2000)
wave_anal = WaveletAnalysis(data, wavelet=wavelets.Ricker())
widths = wave_anal.scales[::-1]
morlet = scipy.signal.morlet
cwt = cwt(data, morlet, widths)
fft_cwt = fft_cwt(data, morlet, widths)
npt.assert_array_almost_equal(cwt, fft_cwt, decimal=13)
def test_random_x(self):
'''
Automatically generate a number of test cases for DFT() and compare results to numpy.fft.ftt
This will generate and test input lists of length 2 to max_N a total of 10 times each
'''
max_N=20
#for x of length 2 to max_N (inclusive)
for N in range(2,max_N+1):
#generate and test x ten times
for t in range(0,10):
#randomly generate x
x = []
for i in range(0,N):
x.append((random()-0.5)*2+(random()-0.5)*2j)
#test DFT by comparing to fft.fft out to 6 decimal places
testing.assert_array_almost_equal(DFT(x),fft.fft(x), err_msg='Your results do not agree with numpy.fft.fft',verbose=True)
#plot_rand(x)
def test_timings(self):
'''Test that timings gives us the right results.'''
timings = analyze.timings(self.samples)
pdt.assert_index_equal(
timings.index,
self.samples.index)
self.assertSetEqual(
set(timings.columns),
set(self.samples.columns.get_level_values(0)))
self.assertEqual(
len(timings.columns),
len(self.samples.columns) / 2)
npt.assert_array_almost_equal(
timings['fn1'],
1.1)
npt.assert_array_almost_equal(
timings['fn2'],
1.5)
def test_isolate(self):
'''Test that isolate gives us the right results.'''
isolated = analyze.isolate(self.samples)
pdt.assert_index_equal(
isolated.index,
self.samples.index)
pdt.assert_index_equal(
isolated.columns,
self.samples.columns)
pdt.assert_frame_equal(
analyze.timings(isolated),
analyze.timings(self.samples))
npt.assert_array_almost_equal(
isolated['fn1']['begin'],
np.arange(20) * 1.1)
npt.assert_array_almost_equal(
isolated['fn1']['end'],
1.1 + (np.arange(20) * 1.1))
npt.assert_array_almost_equal(
isolated['fn2']['begin'],
np.arange(20) * 1.5)
npt.assert_array_almost_equal(
isolated['fn2']['end'],
1.5 + (np.arange(20) * 1.5))
def test_predict_dki():
with nbtmp.InTemporaryDirectory() as tmpdir:
fbval = op.join(tmpdir, 'dki.bval')
fbvec = op.join(tmpdir, 'dki.bvec')
fdata = op.join(tmpdir, 'dki.nii.gz')
make_dki_data(fbval, fbvec, fdata)
cmd1 = ["pyAFQ_dki", "-d", fdata, "-l", fbval, "-c", fbvec,
"-o", tmpdir]
out = runner.run_command(cmd1)
npt.assert_equal(out[0], 0)
# Get expected values
fparams = op.join(tmpdir, "dki_params.nii.gz")
cmd2 = ["pyAFQ_dki_predict", "-p", fparams, "-l", fbval, "-c", fbvec,
"-o", tmpdir, '-b', '0']
out = runner.run_command(cmd2)
npt.assert_equal(out[0], 0)
pred = nib.load(op.join(tmpdir, "dki_prediction.nii.gz")).get_data()
data = nib.load(op.join(tmpdir, "dki.nii.gz")).get_data()
npt.assert_array_almost_equal(pred, data)
def test_predict_dti():
with nbtmp.InTemporaryDirectory() as tmpdir:
fbval = op.join(tmpdir, 'dti.bval')
fbvec = op.join(tmpdir, 'dti.bvec')
fdata = op.join(tmpdir, 'dti.nii.gz')
make_dti_data(fbval, fbvec, fdata)
cmd1 = ["pyAFQ_dti", "-d", fdata, "-l", fbval, "-c", fbvec,
"-o", tmpdir]
out = runner.run_command(cmd1)
npt.assert_equal(out[0], 0)
# Get expected values
fparams = op.join(tmpdir, "dti_params.nii.gz")
cmd2 = ["pyAFQ_dti_predict", "-p", fparams, "-l", fbval, "-c", fbvec,
"-o", tmpdir, '-b', '0']
out = runner.run_command(cmd2)
npt.assert_equal(out[0], 0)
pred = nib.load(op.join(tmpdir, "dti_prediction.nii.gz")).get_data()
data = nib.load(op.join(tmpdir, "dti.nii.gz")).get_data()
npt.assert_array_almost_equal(pred, data)
def test_hessian(self):
"""Tests hessian calculation
"""
indep = np.arange(10)
resp_mat = np.array([(1.02 * 2 * indep)**2 + 1 * indep -
(2 * indep)**2 - indep,
(2 * indep)**2 + 1.02 * indep -
(2 * indep)**2 - indep])
inp_mat = np.array([[0.02 * 2, 0], [0, 0.02]])
init_sens = np.linalg.lstsq(inp_mat, resp_mat)[0].T
resp_mat = np.array([(1.02**2 * 2 * indep)**2 + indep -
(1.02 * 2 * indep)**2 - indep,
(1.02 * 2 * indep)**2 + 1.02 * indep -
(1.02 * 2 * indep)**2 - indep])
inp_mat = np.array([[0.02 * 1.02 * 2, 0], [0, 0.02]])
step1_sens = np.linalg.lstsq(inp_mat, resp_mat)[0].T
resp_mat = np.array([(1.02 * 2 * indep)**2 + 1.02 * indep -
(2 * indep)**2 - 1.02 * indep,
(2 * indep)**2 + 1.02**2 * indep -
(2 * indep)**2 - 1.02 * indep])
inp_mat = np.array([[0.02 * 2, 0], [0, 0.02 * 1.02]])
step2_sens = np.linalg.lstsq(inp_mat, resp_mat)[0].T
true_hess = np.zeros((2, 10, 2))
true_hess[0, :, :] = (step1_sens - init_sens) / (0.02 * 2.0)
true_hess[1, :, :] = (step2_sens - init_sens) / (0.02 * 1.0)
hessian = self.bayes.get_hessian()
npt.assert_array_almost_equal(hessian['simple'], true_hess,
decimal=8,
err_msg='Hessian not as predicted')
def test_mult_simPQ(self):
"""Test the sensitivity calculation for multiple experiments
"""
new = self.bayes.update(simulations={'simple1': [self.sim1, self.exp1],
'simple2': [self.sim1, self.exp1]})
initial_data = new.get_data()
sens_matrix = new._get_sens(initial_data=initial_data)
P, q = new._get_sim_pq(initial_data, sens_matrix)
sigma = self.exp1.get_sigma()
P_true = np.dot(np.dot(sens_matrix['simple1'].T,
inv(sigma)),
sens_matrix['simple1'])
P_true += P_true
npt.assert_array_almost_equal(P, P_true, decimal=8)
epsilon = new.simulations['simple1']['exp'].\
compare(initial_data['simple1'])
q_true = -np.dot(np.dot(epsilon, inv(sigma)), sens_matrix['simple1'])
q_true += q_true
npt.assert_array_almost_equal(q, q_true, decimal=8)
