# MEUL with a non-parametric homeostasis

In this notebook, we will study how homeostasis (cooperation) may be an essential ingredient to this algorithm working on a winner-take-all basis (competition). This extension has been published as Perrinet, Neural Computation (2010) (see https://laurentperrinet.github.io/publication/perrinet-10-shl ). Compared to other posts, such as this previous post, we improve the code to not depend on any parameter (namely the `C`parameter of the rescaling function). For that, we will use a non-parametric approach based on the use of cumulative histograms.

This is joint work with Victor Boutin and Angelo Francisioni. See also the other posts on unsupervised learning.

In [1]:
```import matplotlib.pyplot as plt
import matplotlib
%matplotlib inline
import numpy as np
np.set_printoptions(formatter = dict( float = lambda x: "%.3g" % x ), precision=3, suppress=True, threshold=np.inf)
from shl_scripts.shl_experiments import SHL
import time
```
In [2]:
```DEBUG = True
DEBUG = False
if not DEBUG:
tag = '2017-12-01_Testing_COMPs'
DEBUG_DOWNSCALE = 1
else:
tag = '2017-12-01_Testing_COMPs-DEBUG'
DEBUG_DOWNSCALE = 10

seed = 42
nb_quant = 512
C = 0.
do_sym = False
verbose = False
i_sample = 42

from shl_scripts.shl_experiments import SHL
matname = tag + '_nohomeo'
shl = SHL(DEBUG_DOWNSCALE=DEBUG_DOWNSCALE,
datapath='/tmp/database', seed=seed,
eta=0.05, verbose=2, record_each=50, n_iter=1000,
eta_homeo=0., alpha_homeo=1.,
do_sym=do_sym, nb_quant=nb_quant)
data = shl.get_data(matname=matname)
```
```loading the data called : data_cache/2017-12-01_Testing_COMPs_nohomeo_data
```
In [3]:
```test_size = shl.batch_size # data.shape[0]//2
data_training = data[:test_size, :]
data_test = data[test_size:,:]
if DEBUG:
test_size = data.shape[0]//20
data_training = data[:(data.shape[0]-test_size),:].copy()
data_test = data[:test_size, :].copy()
```
In [4]:
```dico_partial_learning = shl.learn_dico(data=data_training, matname=matname)
```
```loading the dico called : data_cache/2017-12-01_Testing_COMPs_nohomeo_dico.pkl
```

We start off by using a short learning with no homeostasis such that we end up with a unbalanced dictionary:

In [5]:
```fig, ax = shl.show_dico(dico_partial_learning, data=data, title=matname)
fig.show()
fig, ax = shl.time_plot(dico_partial_learning, variable='prob_active');
fig.show()
```

### MP classique¶

In [6]:
```n_samples, n_pixels = data_test.shape
n_dictionary, n_pixels = dico_partial_learning.dictionary.shape
norm_each_filter = np.sqrt(np.sum(dico_partial_learning.dictionary**2, axis=1))
dico_partial_learning.dictionary /= norm_each_filter[:,np.newaxis]

sparse_code_mp = shl.code(data_test, dico_partial_learning, matname=matname)

from shl_scripts.shl_tools import plot_proba_histogram
fig, ax = plot_proba_histogram(sparse_code_mp)
```
```loading the code called : data_cache/2017-12-01_Testing_COMPs_nohomeo_coding.npy
```

### COMP : learning rescaling on sparse coefficients - basic principle¶

In [7]:
```n_samples, n_dictionary = sparse_code_mp.shape
N = n_samples * n_dictionary
q = np.zeros_like(sparse_code_mp)
for i in range(i_sample, i_sample+1): #range(n_samples):
for k in range(n_dictionary):
if sparse_code_mp[i, k] > 0:
q[i, k] = np.sum(sparse_code_mp[i, k]>sparse_code_mp.ravel())/N
if q[i, k]==0:
print(k, i, sparse_code_mp[i, k])
if i < 22 and k < 20:
print(i, k, 'Raw value=', sparse_code_mp[i, k],
' is transformed into=', q[i, k])
```
In [8]:
```def rescaling(code, do_sym=do_sym):
if do_sym:
code = np.abs(code)
else:
code *= code>0
n_samples, n_dictionary = code.shape
N = n_samples * n_dictionary
q = np.zeros_like(code)
for i in range(n_samples):
for k in range(n_dictionary):
if code[i, k] > 0:
q[i, k] = np.sum(code[i, k]>code.ravel())/N
return q
```

Testing on one sample:

In [9]:
```%%time
sparse_code_mp_ = rescaling(sparse_code_mp[i_sample:(i_sample+1), :], do_sym=do_sym)
```
```CPU times: user 275 µs, sys: 10 µs, total: 285 µs
Wall time: 281 µs
```
In [10]:
```def rescaling(code, do_sym=do_sym):
if do_sym: code = np.abs(code)
n_samples, n_dictionary = code.shape
N = n_samples * n_dictionary
q = np.sum(code[:, :, np.newaxis] > code.ravel()[np.newaxis, np.newaxis, :], axis=-1)
return q
```
In [11]:
```%%time
sparse_code_mp_ = rescaling(sparse_code_mp[i_sample:(i_sample+1), :], do_sym=do_sym)
```
```CPU times: user 590 µs, sys: 545 µs, total: 1.13 ms
Wall time: 748 µs
```
In [12]:
```def get_P_cum(code, nb_quant=512, do_rescale=True, do_sym=do_sym):
if do_rescale:
p_c = rescaling(code, do_sym=do_sym)

n_samples, nb_filter = code.shape
code_bins = np.linspace(0., 1, nb_quant, endpoint=True)
P_cum = np.zeros((nb_filter, nb_quant))
for i in range(nb_filter):
p, bins = np.histogram(p_c[:, i], bins=code_bins, density=True)
p /= p.sum()
P_cum[i, :] = np.hstack((0, np.cumsum(p)))
return P_cum
```
In [13]:
```%%time
P_cum = get_P_cum(sparse_code_mp[i_sample:(i_sample+1), :], nb_quant=nb_quant, do_sym=do_sym)
```
```CPU times: user 65.6 ms, sys: 2.27 ms, total: 67.9 ms
Wall time: 68.1 ms
```
P_cum = get_P_cum(sparse_code_mp[i_sample:(i_sample+1), :], nb_quant=nb_quant, do_sym=do_sym) from shl_scripts.shl_tools import plot_P_cum fig, ax = plot_P_cum(P_cum, verbose=False, alpha=.05) ax.plot([0, 1], [0, 1], 'r--') #ax.set_ylim(0.8, 1.01);ax.set_xlim(0.8, 1.01);

### COMP : learning modulations on linear coefficients using ranks¶

Sparse coefficients are distributed according to:

In [14]:
```print('Min=', np.min(sparse_code_mp), 'Max=', np.max(sparse_code_mp), 'Shape=', sparse_code_mp.shape)
fig, ax = plt.subplots(1, 1, figsize=(8, 5))
ax.plot(np.sort(sparse_code_mp.ravel()))
ax.set_xlabel('rank'); ax.set_ylabel('coefficient');
```
```Min= 0.0 Max= 36.576175154 Shape= (81792, 324)
```

and linear (rectified) coefficients to:

In [15]:
```corr = (data_test @ dico_partial_learning.dictionary.T)
corr *= corr>0

print('Min=', np.min(corr), 'Max=', np.max(corr), 'Shape=', corr.shape)

fig, ax = plt.subplots(1, 1, figsize=(8, 5))
ax.plot(np.sort(corr.ravel()))
ax.set_xlabel('rank'); ax.set_ylabel('coefficient');
```
```Min= -0.0 Max= 36.576175154 Shape= (81792, 324)
```

The rescaling functions aims at achieving a better numerical stability, such that the quantification won't affect the learning. But it fundamentally does not change the learning algorithm. The later choice (linear rectified coefficients) seems to be a more robust choice and is similar to the prior choice of an exponential function in a past implementation.

In [16]:
```def get_rescaling(code, nb_quant=nb_quant, do_sym=do_sym, verbose=False):
#     if do_sym:
#         code = np.abs(code)
#     else:
#         code *= code>0
sorted_coeffs = np.sort(code.ravel())
indices = [int(q*(sorted_coeffs.size-1) ) for q in np.linspace(0, 1, nb_quant, endpoint=True)]
C = sorted_coeffs[indices]
if verbose:
print ('At indices (ranks) ', indices)
print ('the coefficients are ', C )
return C

C = get_rescaling(corr, nb_quant=nb_quant, do_sym=do_sym, verbose=True)
```
```At indices (ranks)  [0, 51860, 103720, 155580, 207441, 259301, 311161, 363022, 414882, 466742, 518602, 570463, 622323, 674183, 726044, 777904, 829764, 881624, 933485, 985345, 1037205, 1089066, 1140926, 1192786, 1244646, 1296507, 1348367, 1400227, 1452088, 1503948, 1555808, 1607668, 1659529, 1711389, 1763249, 1815110, 1866970, 1918830, 1970690, 2022551, 2074411, 2126271, 2178132, 2229992, 2281852, 2333712, 2385573, 2437433, 2489293, 2541154, 2593014, 2644874, 2696734, 2748595, 2800455, 2852315, 2904176, 2956036, 3007896, 3059756, 3111617, 3163477, 3215337, 3267198, 3319058, 3370918, 3422778, 3474639, 3526499, 3578359, 3630220, 3682080, 3733940, 3785801, 3837661, 3889521, 3941381, 3993242, 4045102, 4096962, 4148823, 4200683, 4252543, 4304403, 4356264, 4408124, 4459984, 4511845, 4563705, 4615565, 4667425, 4719286, 4771146, 4823006, 4874867, 4926727, 4978587, 5030447, 5082308, 5134168, 5186028, 5237889, 5289749, 5341609, 5393469, 5445330, 5497190, 5549050, 5600911, 5652771, 5704631, 5756491, 5808352, 5860212, 5912072, 5963933, 6015793, 6067653, 6119513, 6171374, 6223234, 6275094, 6326955, 6378815, 6430675, 6482535, 6534396, 6586256, 6638116, 6689977, 6741837, 6793697, 6845557, 6897418, 6949278, 7001138, 7052999, 7104859, 7156719, 7208579, 7260440, 7312300, 7364160, 7416021, 7467881, 7519741, 7571602, 7623462, 7675322, 7727182, 7779043, 7830903, 7882763, 7934624, 7986484, 8038344, 8090204, 8142065, 8193925, 8245785, 8297646, 8349506, 8401366, 8453226, 8505087, 8556947, 8608807, 8660668, 8712528, 8764388, 8816248, 8868109, 8919969, 8971829, 9023690, 9075550, 9127410, 9179270, 9231131, 9282991, 9334851, 9386712, 9438572, 9490432, 9542292, 9594153, 9646013, 9697873, 9749734, 9801594, 9853454, 9905314, 9957175, 10009035, 10060895, 10112756, 10164616, 10216476, 10268336, 10320197, 10372057, 10423917, 10475778, 10527638, 10579498, 10631358, 10683219, 10735079, 10786939, 10838800, 10890660, 10942520, 10994380, 11046241, 11098101, 11149961, 11201822, 11253682, 11305542, 11357403, 11409263, 11461123, 11512983, 11564844, 11616704, 11668564, 11720425, 11772285, 11824145, 11876005, 11927866, 11979726, 12031586, 12083447, 12135307, 12187167, 12239027, 12290888, 12342748, 12394608, 12446469, 12498329, 12550189, 12602049, 12653910, 12705770, 12757630, 12809491, 12861351, 12913211, 12965071, 13016932, 13068792, 13120652, 13172513, 13224373, 13276233, 13328093, 13379954, 13431814, 13483674, 13535535, 13587395, 13639255, 13691115, 13742976, 13794836, 13846696, 13898557, 13950417, 14002277, 14054137, 14105998, 14157858, 14209718, 14261579, 14313439, 14365299, 14417159, 14469020, 14520880, 14572740, 14624601, 14676461, 14728321, 14780181, 14832042, 14883902, 14935762, 14987623, 15039483, 15091343, 15143204, 15195064, 15246924, 15298784, 15350645, 15402505, 15454365, 15506226, 15558086, 15609946, 15661806, 15713667, 15765527, 15817387, 15869248, 15921108, 15972968, 16024828, 16076689, 16128549, 16180409, 16232270, 16284130, 16335990, 16387850, 16439711, 16491571, 16543431, 16595292, 16647152, 16699012, 16750872, 16802733, 16854593, 16906453, 16958314, 17010174, 17062034, 17113894, 17165755, 17217615, 17269475, 17321336, 17373196, 17425056, 17476916, 17528777, 17580637, 17632497, 17684358, 17736218, 17788078, 17839938, 17891799, 17943659, 17995519, 18047380, 18099240, 18151100, 18202960, 18254821, 18306681, 18358541, 18410402, 18462262, 18514122, 18565982, 18617843, 18669703, 18721563, 18773424, 18825284, 18877144, 18929004, 18980865, 19032725, 19084585, 19136446, 19188306, 19240166, 19292027, 19343887, 19395747, 19447607, 19499468, 19551328, 19603188, 19655049, 19706909, 19758769, 19810629, 19862490, 19914350, 19966210, 20018071, 20069931, 20121791, 20173651, 20225512, 20277372, 20329232, 20381093, 20432953, 20484813, 20536673, 20588534, 20640394, 20692254, 20744115, 20795975, 20847835, 20899695, 20951556, 21003416, 21055276, 21107137, 21158997, 21210857, 21262717, 21314578, 21366438, 21418298, 21470159, 21522019, 21573879, 21625739, 21677600, 21729460, 21781320, 21833181, 21885041, 21936901, 21988761, 22040622, 22092482, 22144342, 22196203, 22248063, 22299923, 22351783, 22403644, 22455504, 22507364, 22559225, 22611085, 22662945, 22714806, 22766666, 22818526, 22870386, 22922247, 22974107, 23025967, 23077828, 23129688, 23181548, 23233408, 23285269, 23337129, 23388989, 23440850, 23492710, 23544570, 23596430, 23648291, 23700151, 23752011, 23803872, 23855732, 23907592, 23959452, 24011313, 24063173, 24115033, 24166894, 24218754, 24270614, 24322474, 24374335, 24426195, 24478055, 24529916, 24581776, 24633636, 24685496, 24737357, 24789217, 24841077, 24892938, 24944798, 24996658, 25048518, 25100379, 25152239, 25204099, 25255960, 25307820, 25359680, 25411540, 25463401, 25515261, 25567121, 25618982, 25670842, 25722702, 25774562, 25826423, 25878283, 25930143, 25982004, 26033864, 26085724, 26137584, 26189445, 26241305, 26293165, 26345026, 26396886, 26448746, 26500607]
the coefficients are  [-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 0.000709 0.00149 0.00227 0.00307 0.00386 0.00461
0.00537 0.00613 0.00693 0.00772 0.00851 0.00929 0.0101 0.0109 0.0117
0.0126 0.0134 0.0142 0.0151 0.016 0.0169 0.0178 0.0187 0.0196 0.0206
0.0215 0.0225 0.0235 0.0245 0.0255 0.0265 0.0275 0.0286 0.0297 0.0308
0.0321 0.0333 0.0345 0.0358 0.0371 0.0384 0.0398 0.0413 0.0428 0.0443
0.046 0.0477 0.0494 0.0511 0.0529 0.0547 0.0567 0.0587 0.0607 0.0628 0.065
0.0671 0.0693 0.0717 0.0741 0.0766 0.0792 0.0819 0.0846 0.0873 0.0901
0.093 0.0959 0.0989 0.102 0.105 0.108 0.111 0.115 0.118 0.121 0.125 0.128
0.132 0.135 0.139 0.143 0.146 0.15 0.154 0.158 0.162 0.166 0.17 0.174
0.178 0.183 0.187 0.191 0.195 0.2 0.204 0.209 0.214 0.218 0.223 0.228
0.233 0.238 0.243 0.248 0.253 0.258 0.263 0.269 0.274 0.28 0.285 0.291
0.296 0.302 0.308 0.314 0.32 0.326 0.332 0.338 0.344 0.351 0.357 0.363
0.37 0.377 0.384 0.39 0.397 0.405 0.412 0.419 0.426 0.434 0.441 0.449
0.457 0.464 0.472 0.48 0.489 0.497 0.505 0.514 0.522 0.531 0.54 0.549
0.558 0.567 0.577 0.586 0.596 0.606 0.615 0.626 0.636 0.646 0.657 0.667
0.678 0.689 0.701 0.712 0.723 0.735 0.747 0.759 0.772 0.784 0.797 0.81
0.823 0.836 0.85 0.864 0.878 0.893 0.907 0.922 0.937 0.952 0.968 0.984 1
1.02 1.03 1.05 1.07 1.09 1.11 1.12 1.14 1.16 1.18 1.2 1.22 1.25 1.27 1.29
1.31 1.34 1.36 1.39 1.41 1.44 1.46 1.49 1.52 1.55 1.58 1.61 1.64 1.67 1.7
1.74 1.77 1.81 1.85 1.89 1.93 1.97 2.01 2.06 2.1 2.15 2.2 2.26 2.31 2.37
2.43 2.5 2.57 2.64 2.71 2.8 2.88 2.98 3.08 3.19 3.31 3.44 3.59 3.75 3.94
4.16 4.42 4.74 5.17 5.8 6.91 36.6]
```
In [17]:
```def rescaling(code, C, do_sym=do_sym):
#     if do_sym:
#         code = np.abs(code)
#     else:
#         code *= code>0
n_samples, n_dictionary = code.shape
N = n_samples * n_dictionary
q = np.zeros_like(code)
code_bins = np.linspace(0., 1, C.size, endpoint=True)

for i in range(n_samples):
for k in range(n_dictionary):
if code[i, k] > 0:
q[i, k] = np.interp(code[i, k], C, code_bins)
return q

rescaled_code_ = rescaling(corr[i_sample:(i_sample+1), :], C=C, do_sym=do_sym)
```
In [18]:
```%%time
rescaled_code_ = rescaling(corr[i_sample:(i_sample+1), :], C=C, do_sym=do_sym)
```
```CPU times: user 1.89 ms, sys: 1.36 ms, total: 3.25 ms
Wall time: 1.98 ms
```

This can be vectorized:

In [19]:
```def rescaling(code, C, do_sym=do_sym):
#     if do_sym:
#         code = np.abs(code)
#     else:
#         code *= code>0
code_bins = np.linspace(0., 1, C.size, endpoint=True)
return np.interp(code, C, code_bins) * (code > 0.)

print('on ', corr.shape[0], ' samples')
rescaled_code_ = rescaling(corr, C=C, do_sym=do_sym)
```
```on  81792  samples
```
In [20]:
```%%time
sparse_code_mp_ = rescaling(sparse_code_mp, C=C, do_sym=do_sym)
```
```CPU times: user 444 ms, sys: 117 ms, total: 561 ms
Wall time: 635 ms
```

We notice also that we use half of the vector for negative coefficients. This could be improved as:

In [21]:
```def rescaling(code, C, do_sym=do_sym):
#     if do_sym:
#         code = np.abs(code)
#     else:
#         code *= code>0
p_c = np.zeros_like(code)
ind_nz = code>0.
code_bins = np.linspace(0., 1, C.size, endpoint=True)
p_c[ind_nz] = np.interp(code[ind_nz], C, code_bins)
return p_c # * (code > 0.)

rescaled_code__ = rescaling(corr, C=C, do_sym=do_sym)
assert(np.array_equal(rescaled_code_, rescaled_code__))
```
In [22]:
```%%time
sparse_code_mp_ = rescaling(sparse_code_mp, C=C, do_sym=do_sym)
```
```CPU times: user 246 ms, sys: 131 ms, total: 377 ms
Wall time: 492 ms
```

Finally, one can get the `P_cum` look-up-tables:

In [23]:
```def get_P_cum(code, nb_quant, C=C, do_sym=do_sym):
p_c = rescaling(code, C=C, do_sym=do_sym)

n_samples, nb_filter = code.shape
code_bins = np.linspace(0., 1, nb_quant, endpoint=True)
P_cum = np.zeros((nb_filter, nb_quant))
for i in range(nb_filter):
p, bins = np.histogram(p_c[:, i], bins=code_bins, density=True)
p /= p.sum()
P_cum[i, :] = np.hstack((0., np.cumsum(p)))
return P_cum
```

computing histograms is rather fast:

In [24]:
```%%time
P_cum = get_P_cum(sparse_code_mp, nb_quant=nb_quant, C=C)
```
```CPU times: user 1.23 s, sys: 126 ms, total: 1.36 s
Wall time: 1.44 s
```
In [25]:
```C = get_rescaling(corr, nb_quant=nb_quant, do_sym=do_sym, verbose=True)
C_vec_init = C.copy()
fig=plt.figure(figsize=(13, 8))
ax = plt.subplot(111)
x = np.linspace(0, 30, 100)
ax.plot(x, np.ones_like(x), '--')
ax.plot(x, rescaling(x, C=C, do_sym=do_sym))
ax.set_xlabel('raw coefficient')
ax.set_ylabel('rescaled coefficient')
ax.set_xscale('log');
```
```At indices (ranks)  [0, 51860, 103720, 155580, 207441, 259301, 311161, 363022, 414882, 466742, 518602, 570463, 622323, 674183, 726044, 777904, 829764, 881624, 933485, 985345, 1037205, 1089066, 1140926, 1192786, 1244646, 1296507, 1348367, 1400227, 1452088, 1503948, 1555808, 1607668, 1659529, 1711389, 1763249, 1815110, 1866970, 1918830, 1970690, 2022551, 2074411, 2126271, 2178132, 2229992, 2281852, 2333712, 2385573, 2437433, 2489293, 2541154, 2593014, 2644874, 2696734, 2748595, 2800455, 2852315, 2904176, 2956036, 3007896, 3059756, 3111617, 3163477, 3215337, 3267198, 3319058, 3370918, 3422778, 3474639, 3526499, 3578359, 3630220, 3682080, 3733940, 3785801, 3837661, 3889521, 3941381, 3993242, 4045102, 4096962, 4148823, 4200683, 4252543, 4304403, 4356264, 4408124, 4459984, 4511845, 4563705, 4615565, 4667425, 4719286, 4771146, 4823006, 4874867, 4926727, 4978587, 5030447, 5082308, 5134168, 5186028, 5237889, 5289749, 5341609, 5393469, 5445330, 5497190, 5549050, 5600911, 5652771, 5704631, 5756491, 5808352, 5860212, 5912072, 5963933, 6015793, 6067653, 6119513, 6171374, 6223234, 6275094, 6326955, 6378815, 6430675, 6482535, 6534396, 6586256, 6638116, 6689977, 6741837, 6793697, 6845557, 6897418, 6949278, 7001138, 7052999, 7104859, 7156719, 7208579, 7260440, 7312300, 7364160, 7416021, 7467881, 7519741, 7571602, 7623462, 7675322, 7727182, 7779043, 7830903, 7882763, 7934624, 7986484, 8038344, 8090204, 8142065, 8193925, 8245785, 8297646, 8349506, 8401366, 8453226, 8505087, 8556947, 8608807, 8660668, 8712528, 8764388, 8816248, 8868109, 8919969, 8971829, 9023690, 9075550, 9127410, 9179270, 9231131, 9282991, 9334851, 9386712, 9438572, 9490432, 9542292, 9594153, 9646013, 9697873, 9749734, 9801594, 9853454, 9905314, 9957175, 10009035, 10060895, 10112756, 10164616, 10216476, 10268336, 10320197, 10372057, 10423917, 10475778, 10527638, 10579498, 10631358, 10683219, 10735079, 10786939, 10838800, 10890660, 10942520, 10994380, 11046241, 11098101, 11149961, 11201822, 11253682, 11305542, 11357403, 11409263, 11461123, 11512983, 11564844, 11616704, 11668564, 11720425, 11772285, 11824145, 11876005, 11927866, 11979726, 12031586, 12083447, 12135307, 12187167, 12239027, 12290888, 12342748, 12394608, 12446469, 12498329, 12550189, 12602049, 12653910, 12705770, 12757630, 12809491, 12861351, 12913211, 12965071, 13016932, 13068792, 13120652, 13172513, 13224373, 13276233, 13328093, 13379954, 13431814, 13483674, 13535535, 13587395, 13639255, 13691115, 13742976, 13794836, 13846696, 13898557, 13950417, 14002277, 14054137, 14105998, 14157858, 14209718, 14261579, 14313439, 14365299, 14417159, 14469020, 14520880, 14572740, 14624601, 14676461, 14728321, 14780181, 14832042, 14883902, 14935762, 14987623, 15039483, 15091343, 15143204, 15195064, 15246924, 15298784, 15350645, 15402505, 15454365, 15506226, 15558086, 15609946, 15661806, 15713667, 15765527, 15817387, 15869248, 15921108, 15972968, 16024828, 16076689, 16128549, 16180409, 16232270, 16284130, 16335990, 16387850, 16439711, 16491571, 16543431, 16595292, 16647152, 16699012, 16750872, 16802733, 16854593, 16906453, 16958314, 17010174, 17062034, 17113894, 17165755, 17217615, 17269475, 17321336, 17373196, 17425056, 17476916, 17528777, 17580637, 17632497, 17684358, 17736218, 17788078, 17839938, 17891799, 17943659, 17995519, 18047380, 18099240, 18151100, 18202960, 18254821, 18306681, 18358541, 18410402, 18462262, 18514122, 18565982, 18617843, 18669703, 18721563, 18773424, 18825284, 18877144, 18929004, 18980865, 19032725, 19084585, 19136446, 19188306, 19240166, 19292027, 19343887, 19395747, 19447607, 19499468, 19551328, 19603188, 19655049, 19706909, 19758769, 19810629, 19862490, 19914350, 19966210, 20018071, 20069931, 20121791, 20173651, 20225512, 20277372, 20329232, 20381093, 20432953, 20484813, 20536673, 20588534, 20640394, 20692254, 20744115, 20795975, 20847835, 20899695, 20951556, 21003416, 21055276, 21107137, 21158997, 21210857, 21262717, 21314578, 21366438, 21418298, 21470159, 21522019, 21573879, 21625739, 21677600, 21729460, 21781320, 21833181, 21885041, 21936901, 21988761, 22040622, 22092482, 22144342, 22196203, 22248063, 22299923, 22351783, 22403644, 22455504, 22507364, 22559225, 22611085, 22662945, 22714806, 22766666, 22818526, 22870386, 22922247, 22974107, 23025967, 23077828, 23129688, 23181548, 23233408, 23285269, 23337129, 23388989, 23440850, 23492710, 23544570, 23596430, 23648291, 23700151, 23752011, 23803872, 23855732, 23907592, 23959452, 24011313, 24063173, 24115033, 24166894, 24218754, 24270614, 24322474, 24374335, 24426195, 24478055, 24529916, 24581776, 24633636, 24685496, 24737357, 24789217, 24841077, 24892938, 24944798, 24996658, 25048518, 25100379, 25152239, 25204099, 25255960, 25307820, 25359680, 25411540, 25463401, 25515261, 25567121, 25618982, 25670842, 25722702, 25774562, 25826423, 25878283, 25930143, 25982004, 26033864, 26085724, 26137584, 26189445, 26241305, 26293165, 26345026, 26396886, 26448746, 26500607]
the coefficients are  [-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0 -0
-0 -0 -0 -0 -0 -0 -0 -0 0.000709 0.00149 0.00227 0.00307 0.00386 0.00461
0.00537 0.00613 0.00693 0.00772 0.00851 0.00929 0.0101 0.0109 0.0117
0.0126 0.0134 0.0142 0.0151 0.016 0.0169 0.0178 0.0187 0.0196 0.0206
0.0215 0.0225 0.0235 0.0245 0.0255 0.0265 0.0275 0.0286 0.0297 0.0308
0.0321 0.0333 0.0345 0.0358 0.0371 0.0384 0.0398 0.0413 0.0428 0.0443
0.046 0.0477 0.0494 0.0511 0.0529 0.0547 0.0567 0.0587 0.0607 0.0628 0.065
0.0671 0.0693 0.0717 0.0741 0.0766 0.0792 0.0819 0.0846 0.0873 0.0901
0.093 0.0959 0.0989 0.102 0.105 0.108 0.111 0.115 0.118 0.121 0.125 0.128
0.132 0.135 0.139 0.143 0.146 0.15 0.154 0.158 0.162 0.166 0.17 0.174
0.178 0.183 0.187 0.191 0.195 0.2 0.204 0.209 0.214 0.218 0.223 0.228
0.233 0.238 0.243 0.248 0.253 0.258 0.263 0.269 0.274 0.28 0.285 0.291
0.296 0.302 0.308 0.314 0.32 0.326 0.332 0.338 0.344 0.351 0.357 0.363
0.37 0.377 0.384 0.39 0.397 0.405 0.412 0.419 0.426 0.434 0.441 0.449
0.457 0.464 0.472 0.48 0.489 0.497 0.505 0.514 0.522 0.531 0.54 0.549
0.558 0.567 0.577 0.586 0.596 0.606 0.615 0.626 0.636 0.646 0.657 0.667
0.678 0.689 0.701 0.712 0.723 0.735 0.747 0.759 0.772 0.784 0.797 0.81
0.823 0.836 0.85 0.864 0.878 0.893 0.907 0.922 0.937 0.952 0.968 0.984 1
1.02 1.03 1.05 1.07 1.09 1.11 1.12 1.14 1.16 1.18 1.2 1.22 1.25 1.27 1.29
1.31 1.34 1.36 1.39 1.41 1.44 1.46 1.49 1.52 1.55 1.58 1.61 1.64 1.67 1.7
1.74 1.77 1.81 1.85 1.89 1.93 1.97 2.01 2.06 2.1 2.15 2.2 2.26 2.31 2.37
2.43 2.5 2.57 2.64 2.71 2.8 2.88 2.98 3.08 3.19 3.31 3.44 3.59 3.75 3.94
4.16 4.42 4.74 5.17 5.8 6.91 36.6]
```
In [26]:
```P_cum = get_P_cum(rescaling(sparse_code_mp, C=C, do_sym=do_sym), nb_quant=nb_quant, C=C, do_sym=do_sym)
from shl_scripts.shl_tools import plot_P_cum
fig, ax = plot_P_cum(P_cum, verbose=False, alpha=.05)
ax.plot([0, 1], [0, 1], 'r--');
ax.set_ylim(0.85, 1.01);ax.set_xlim(0., 1.);
```

Finally, we notice that the `C` vector is of the same length as each `P_cum` vector and we will integrate it to this array to ease message passing.

### COMP : optimizing the quantile function¶

These `quantile` and `rescaling` functions are implanted in the `shl_scripts`:

In [27]:
```from shl_scripts.shl_encode import quantile, rescaling

corr = (data_test @ dico_partial_learning.dictionary.T)
#corr *= corr>0
print('correlation=', corr[i_sample, :])
print('transformed correlation=', rescaling(corr, C=C, do_sym=do_sym)[i_sample, :])
```
```correlation= [-0.606 0.428 -0.566 -0.457 0.304 0.687 0.55 -1.84 0.43 0.755 -0.286 0.173
0.365 -0.185 0.615 0.263 0.614 0.819 -1.68 1.78 0.0753 -0.536 0.404 -1.77
1.14 0.117 0.625 -0.296 -0.0527 0.591 0.175 -0.639 0.812 0.0353 0.825 1.07
1.97 -0.282 0.187 -0.305 -1.18 0.818 -0.2 1.31 -0.218 -1.13 0.809 0.509
0.069 -0.673 -0.131 -0.165 -0.0833 -1.9 -0.244 -0.414 -0.0485 0.718 -0.154
0.374 0.0455 1.17 -0.048 1.44 0.602 0.00801 -1.67 1.16 0.307 -0.152 1.36
-0.173 0.444 -0.175 0.561 1.35 -0.738 1.36 2.26 0.855 0.666 0.527 -0.456
1.3 0.556 -0.217 -0.28 1.52 0.744 -0.307 -1.02 0.891 6.32e-05 1.09 -1.94
0.568 -0.473 1.05 0.245 0.163 1.13 -0.0183 0.686 0.543 0.379 0.6 1.38
0.456 2.21 -1.09 0.411 -0.872 2.18 -0.301 0.681 1.35 0.568 -0.337 1.66
-1.09 -0.486 -0.598 -0.241 1.52 -0.126 0.42 -0.0964 -1.19 1.09 0.0617
0.542 -0.0428 -1.21 0.549 -0.385 0.2 0.0937 -0.74 -0.845 -0.783 0.0349
0.973 1.13 -0.306 -1.02 0.0624 0.0613 -0.278 0.967 -0.658 0.579 -0.0586
-0.463 -0.726 0.412 -0.0301 0.929 0.765 -0.3 -1.11 -0.222 -0.213 -0.21
0.231 -0.907 -0.117 0.0441 0.918 -0.195 0.603 -0.661 0.419 -0.522 -0.182
1.41 0.414 -0.331 -0.256 0.636 0.407 0.156 0.0501 1.68 0.811 1.4 -0.833
1.01 1.8 0.263 1.28 -0.559 -0.0339 1.19 0.828 -0.212 2.02 0.708 -0.0837
-0.791 -1 0.226 -0.353 1.26 1.18 0.64 0.21 0.61 -0.0966 -1.19 0.58 0.0121
0.283 -1.74 0.911 0.98 -0.381 0.11 0.359 0.502 -0.574 -0.475 0.0498 0.0972
0.775 -1.33 -0.757 -0.0152 1.11 -0.464 0.0499 -1.77 1.81 1.26 0.227 -0.858
0.734 -0.228 -0.602 1.15 1.47 1.52 1.4 1.02 -0.039 0.217 0.654 2.4 0.481
0.161 -0.47 0.626 0.356 0.708 0.68 0.992 -0.364 -0.106 0.948 -0.418 0.712
0.529 -0.207 1.03 -0.315 -0.741 0.481 -0.457 -1.01 -0.638 0.423 -0.278
1.09 1.46 1.95 -0.22 0.571 -0.524 -0.645 -1.01 1.1 0.0848 1.58 0.34 0.143
0.187 0.363 1.24 0.583 -0.577 0.118 -0.611 1.7 -2.21 0.435 0.0115 0.717
0.984 0.666 -1.57 1.36 -0.272 -0.204 -0.0694 0.885 0.157 0.598 0.687
-0.432 0.737 -0.863 0.506 -0.223 -1.04 -0.0878 0.246 -0.237 0.823 -0.371
-0.939 1.31 0.642 -0.593 0.0982 1.84]
transformed correlation= [0 0.768 0 0 0.73 0.823 0.797 0 0.768 0.835 0 0.679 0.75 0 0.81 0.716 0.81
0.845 0 0.932 0.621 0 0.761 0 0.884 0.649 0.812 0 0 0.805 0.68 0 0.844
0.578 0.846 0.876 0.941 0 0.685 0 0 0.845 0 0.9 0 0 0.843 0.788 0.616 0 0
0 0 0 0 0 0 0.829 0 0.753 0.592 0.887 0 0.91 0.808 0.523 0 0.887 0.732 0
0.904 0 0.772 0 0.799 0.904 0 0.904 0.953 0.85 0.82 0.792 0 0.899 0.798 0
0 0.916 0.833 0 0 0.855 0.503 0.879 0 0.801 0 0.874 0.709 0.674 0.883 0
0.823 0.795 0.754 0.807 0.906 0.775 0.951 0 0.763 0 0.95 0 0.822 0.903
0.801 0 0.925 0 0 0 0 0.916 0 0.765 0 0 0.879 0.61 0.795 0 0 0.796 0 0.691
0.635 0 0 0 0.578 0.866 0.883 0 0 0.61 0.609 0 0.865 0 0.803 0 0 0 0.763 0
0.86 0.837 0 0 0 0 0 0.704 0 0 0.591 0.859 0 0.808 0 0.765 0 0 0.908 0.764
0 0 0.814 0.762 0.67 0.598 0.926 0.844 0.907 0 0.87 0.933 0.716 0.897 0 0
0.889 0.846 0 0.944 0.827 0 0 0 0.702 0 0.895 0.889 0.815 0.695 0.809 0 0
0.803 0.533 0.723 0 0.858 0.866 0 0.645 0.748 0.786 0 0 0.597 0.637 0.838
0 0 0 0.881 0 0.597 0 0.934 0.895 0.702 0 0.832 0 0 0.885 0.912 0.916
0.908 0.871 0 0.698 0.817 0.958 0.781 0.673 0 0.812 0.747 0.827 0.822
0.868 0 0 0.862 0 0.828 0.792 0 0.872 0 0 0.781 0 0 0 0.766 0 0.879 0.912
0.94 0 0.801 0 0 0 0.88 0.628 0.92 0.742 0.664 0.685 0.749 0.894 0.804 0
0.649 0 0.927 0 0.769 0.532 0.829 0.867 0.82 0 0.904 0 0 0 0.854 0.671
0.807 0.824 0 0.832 0 0.787 0 0 0 0.71 0 0.845 0 0 0.9 0.815 0 0.638 0.935]
```
In [28]:
```C = get_rescaling(corr, nb_quant=nb_quant, do_sym=do_sym, verbose=True)
fig=plt.figure(figsize=(13, 8))
ax = plt.subplot(111)
x = np.linspace(0, 30, 100)
ax.plot(x, np.ones_like(x), '--')
ax.plot(x, rescaling(x, C=C_vec_init, do_sym=do_sym))
ax.plot(x, rescaling(x, C=.5, do_sym=do_sym))
ax.set_xlabel('raw coefficient')
ax.set_ylabel('rescaled coefficient')
ax.set_xscale('log')
ax.axis('tight');
```
```At indices (ranks)  [0, 51860, 103720, 155580, 207441, 259301, 311161, 363022, 414882, 466742, 518602, 570463, 622323, 674183, 726044, 777904, 829764, 881624, 933485, 985345, 1037205, 1089066, 1140926, 1192786, 1244646, 1296507, 1348367, 1400227, 1452088, 1503948, 1555808, 1607668, 1659529, 1711389, 1763249, 1815110, 1866970, 1918830, 1970690, 2022551, 2074411, 2126271, 2178132, 2229992, 2281852, 2333712, 2385573, 2437433, 2489293, 2541154, 2593014, 2644874, 2696734, 2748595, 2800455, 2852315, 2904176, 2956036, 3007896, 3059756, 3111617, 3163477, 3215337, 3267198, 3319058, 3370918, 3422778, 3474639, 3526499, 3578359, 3630220, 3682080, 3733940, 3785801, 3837661, 3889521, 3941381, 3993242, 4045102, 4096962, 4148823, 4200683, 4252543, 4304403, 4356264, 4408124, 4459984, 4511845, 4563705, 4615565, 4667425, 4719286, 4771146, 4823006, 4874867, 4926727, 4978587, 5030447, 5082308, 5134168, 5186028, 5237889, 5289749, 5341609, 5393469, 5445330, 5497190, 5549050, 5600911, 5652771, 5704631, 5756491, 5808352, 5860212, 5912072, 5963933, 6015793, 6067653, 6119513, 6171374, 6223234, 6275094, 6326955, 6378815, 6430675, 6482535, 6534396, 6586256, 6638116, 6689977, 6741837, 6793697, 6845557, 6897418, 6949278, 7001138, 7052999, 7104859, 7156719, 7208579, 7260440, 7312300, 7364160, 7416021, 7467881, 7519741, 7571602, 7623462, 7675322, 7727182, 7779043, 7830903, 7882763, 7934624, 7986484, 8038344, 8090204, 8142065, 8193925, 8245785, 8297646, 8349506, 8401366, 8453226, 8505087, 8556947, 8608807, 8660668, 8712528, 8764388, 8816248, 8868109, 8919969, 8971829, 9023690, 9075550, 9127410, 9179270, 9231131, 9282991, 9334851, 9386712, 9438572, 9490432, 9542292, 9594153, 9646013, 9697873, 9749734, 9801594, 9853454, 9905314, 9957175, 10009035, 10060895, 10112756, 10164616, 10216476, 10268336, 10320197, 10372057, 10423917, 10475778, 10527638, 10579498, 10631358, 10683219, 10735079, 10786939, 10838800, 10890660, 10942520, 10994380, 11046241, 11098101, 11149961, 11201822, 11253682, 11305542, 11357403, 11409263, 11461123, 11512983, 11564844, 11616704, 11668564, 11720425, 11772285, 11824145, 11876005, 11927866, 11979726, 12031586, 12083447, 12135307, 12187167, 12239027, 12290888, 12342748, 12394608, 12446469, 12498329, 12550189, 12602049, 12653910, 12705770, 12757630, 12809491, 12861351, 12913211, 12965071, 13016932, 13068792, 13120652, 13172513, 13224373, 13276233, 13328093, 13379954, 13431814, 13483674, 13535535, 13587395, 13639255, 13691115, 13742976, 13794836, 13846696, 13898557, 13950417, 14002277, 14054137, 14105998, 14157858, 14209718, 14261579, 14313439, 14365299, 14417159, 14469020, 14520880, 14572740, 14624601, 14676461, 14728321, 14780181, 14832042, 14883902, 14935762, 14987623, 15039483, 15091343, 15143204, 15195064, 15246924, 15298784, 15350645, 15402505, 15454365, 15506226, 15558086, 15609946, 15661806, 15713667, 15765527, 15817387, 15869248, 15921108, 15972968, 16024828, 16076689, 16128549, 16180409, 16232270, 16284130, 16335990, 16387850, 16439711, 16491571, 16543431, 16595292, 16647152, 16699012, 16750872, 16802733, 16854593, 16906453, 16958314, 17010174, 17062034, 17113894, 17165755, 17217615, 17269475, 17321336, 17373196, 17425056, 17476916, 17528777, 17580637, 17632497, 17684358, 17736218, 17788078, 17839938, 17891799, 17943659, 17995519, 18047380, 18099240, 18151100, 18202960, 18254821, 18306681, 18358541, 18410402, 18462262, 18514122, 18565982, 18617843, 18669703, 18721563, 18773424, 18825284, 18877144, 18929004, 18980865, 19032725, 19084585, 19136446, 19188306, 19240166, 19292027, 19343887, 19395747, 19447607, 19499468, 19551328, 19603188, 19655049, 19706909, 19758769, 19810629, 19862490, 19914350, 19966210, 20018071, 20069931, 20121791, 20173651, 20225512, 20277372, 20329232, 20381093, 20432953, 20484813, 20536673, 20588534, 20640394, 20692254, 20744115, 20795975, 20847835, 20899695, 20951556, 21003416, 21055276, 21107137, 21158997, 21210857, 21262717, 21314578, 21366438, 21418298, 21470159, 21522019, 21573879, 21625739, 21677600, 21729460, 21781320, 21833181, 21885041, 21936901, 21988761, 22040622, 22092482, 22144342, 22196203, 22248063, 22299923, 22351783, 22403644, 22455504, 22507364, 22559225, 22611085, 22662945, 22714806, 22766666, 22818526, 22870386, 22922247, 22974107, 23025967, 23077828, 23129688, 23181548, 23233408, 23285269, 23337129, 23388989, 23440850, 23492710, 23544570, 23596430, 23648291, 23700151, 23752011, 23803872, 23855732, 23907592, 23959452, 24011313, 24063173, 24115033, 24166894, 24218754, 24270614, 24322474, 24374335, 24426195, 24478055, 24529916, 24581776, 24633636, 24685496, 24737357, 24789217, 24841077, 24892938, 24944798, 24996658, 25048518, 25100379, 25152239, 25204099, 25255960, 25307820, 25359680, 25411540, 25463401, 25515261, 25567121, 25618982, 25670842, 25722702, 25774562, 25826423, 25878283, 25930143, 25982004, 26033864, 26085724, 26137584, 26189445, 26241305, 26293165, 26345026, 26396886, 26448746, 26500607]
the coefficients are  [-36.2 -6.93 -5.81 -5.19 -4.77 -4.44 -4.17 -3.95 -3.76 -3.59 -3.45 -3.31
-3.19 -3.08 -2.98 -2.89 -2.8 -2.72 -2.64 -2.57 -2.5 -2.43 -2.37 -2.31
-2.26 -2.2 -2.15 -2.1 -2.05 -2.01 -1.97 -1.92 -1.88 -1.84 -1.81 -1.77
-1.73 -1.7 -1.67 -1.63 -1.6 -1.57 -1.54 -1.51 -1.49 -1.46 -1.43 -1.41
-1.38 -1.36 -1.33 -1.31 -1.29 -1.26 -1.24 -1.22 -1.2 -1.18 -1.16 -1.14
-1.12 -1.1 -1.08 -1.07 -1.05 -1.03 -1.01 -0.997 -0.98 -0.964 -0.949 -0.933
-0.918 -0.903 -0.888 -0.874 -0.86 -0.846 -0.832 -0.819 -0.806 -0.793 -0.78
-0.768 -0.755 -0.743 -0.731 -0.72 -0.708 -0.697 -0.686 -0.675 -0.664
-0.653 -0.643 -0.632 -0.622 -0.612 -0.602 -0.592 -0.583 -0.573 -0.564
-0.555 -0.546 -0.537 -0.528 -0.519 -0.511 -0.502 -0.494 -0.486 -0.478
-0.47 -0.462 -0.454 -0.447 -0.439 -0.432 -0.424 -0.417 -0.41 -0.403 -0.396
-0.389 -0.382 -0.375 -0.369 -0.362 -0.356 -0.349 -0.343 -0.337 -0.331
-0.325 -0.319 -0.313 -0.307 -0.301 -0.295 -0.29 -0.284 -0.279 -0.273
-0.268 -0.263 -0.258 -0.252 -0.247 -0.242 -0.238 -0.233 -0.228 -0.223
-0.218 -0.214 -0.209 -0.205 -0.2 -0.196 -0.191 -0.187 -0.183 -0.179 -0.175
-0.17 -0.167 -0.163 -0.159 -0.155 -0.151 -0.147 -0.144 -0.14 -0.136 -0.133
-0.129 -0.126 -0.123 -0.119 -0.116 -0.113 -0.11 -0.107 -0.103 -0.1 -0.0975
-0.0946 -0.0918 -0.0891 -0.0864 -0.0837 -0.0811 -0.0785 -0.076 -0.0736
-0.0712 -0.0689 -0.0667 -0.0645 -0.0623 -0.0603 -0.0583 -0.0564 -0.0546
-0.0528 -0.0511 -0.0495 -0.0479 -0.0463 -0.0447 -0.0433 -0.0419 -0.0405
-0.0391 -0.0378 -0.0366 -0.0353 -0.0341 -0.033 -0.0318 -0.0307 -0.0296
-0.0286 -0.0276 -0.0265 -0.0256 -0.0246 -0.0236 -0.0226 -0.0217 -0.0208
-0.0199 -0.019 -0.0181 -0.0173 -0.0164 -0.0156 -0.0148 -0.014 -0.0132
-0.0124 -0.0116 -0.0108 -0.01 -0.00922 -0.00844 -0.00766 -0.00685 -0.00608
-0.00532 -0.00456 -0.00381 -0.00307 -0.00233 -0.00158 -0.000822 -7.28e-05
0.000709 0.00149 0.00227 0.00307 0.00386 0.00461 0.00537 0.00613 0.00693
0.00772 0.00851 0.00929 0.0101 0.0109 0.0117 0.0126 0.0134 0.0142 0.0151
0.016 0.0169 0.0178 0.0187 0.0196 0.0206 0.0215 0.0225 0.0235 0.0245
0.0255 0.0265 0.0275 0.0286 0.0297 0.0308 0.0321 0.0333 0.0345 0.0358
0.0371 0.0384 0.0398 0.0413 0.0428 0.0443 0.046 0.0477 0.0494 0.0511
0.0529 0.0547 0.0567 0.0587 0.0607 0.0628 0.065 0.0671 0.0693 0.0717
0.0741 0.0766 0.0792 0.0819 0.0846 0.0873 0.0901 0.093 0.0959 0.0989 0.102
0.105 0.108 0.111 0.115 0.118 0.121 0.125 0.128 0.132 0.135 0.139 0.143
0.146 0.15 0.154 0.158 0.162 0.166 0.17 0.174 0.178 0.183 0.187 0.191
0.195 0.2 0.204 0.209 0.214 0.218 0.223 0.228 0.233 0.238 0.243 0.248
0.253 0.258 0.263 0.269 0.274 0.28 0.285 0.291 0.296 0.302 0.308 0.314
0.32 0.326 0.332 0.338 0.344 0.351 0.357 0.363 0.37 0.377 0.384 0.39 0.397
0.405 0.412 0.419 0.426 0.434 0.441 0.449 0.457 0.464 0.472 0.48 0.489
0.497 0.505 0.514 0.522 0.531 0.54 0.549 0.558 0.567 0.577 0.586 0.596
0.606 0.615 0.626 0.636 0.646 0.657 0.667 0.678 0.689 0.701 0.712 0.723
0.735 0.747 0.759 0.772 0.784 0.797 0.81 0.823 0.836 0.85 0.864 0.878
0.893 0.907 0.922 0.937 0.952 0.968 0.984 1 1.02 1.03 1.05 1.07 1.09 1.11
1.12 1.14 1.16 1.18 1.2 1.22 1.25 1.27 1.29 1.31 1.34 1.36 1.39 1.41 1.44
1.46 1.49 1.52 1.55 1.58 1.61 1.64 1.67 1.7 1.74 1.77 1.81 1.85 1.89 1.93
1.97 2.01 2.06 2.1 2.15 2.2 2.26 2.31 2.37 2.43 2.5 2.57 2.64 2.71 2.8
2.88 2.98 3.08 3.19 3.31 3.44 3.59 3.75 3.94 4.16 4.42 4.74 5.17 5.8 6.91
36.6]
```

A sanity check with extremal values:

In [29]:
```stick = np.arange(shl.n_dictionary)*nb_quant
print('Value for ones = ', rescaling(np.inf*np.ones(shl.n_dictionary), C=C))
print('Value for zeros = ', rescaling(np.zeros(shl.n_dictionary), C=C))
```
```Value for ones =  [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]
Value for zeros =  [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
```
In [30]:
```print('Value for ones = ', quantile(P_cum, rescaling(np.inf*np.ones(shl.n_dictionary), C=C), stick))
print('Value for zeros = ', quantile(P_cum, rescaling(np.zeros(shl.n_dictionary), C=C), stick))
```
```Value for ones =  [1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1]
Value for zeros =  [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
```

a `P_cum` array that should not change anything to the way we chose filters:

In [31]:
```P_cum_zeroeffect = np.linspace(0, 1, nb_quant, endpoint=True)[np.newaxis, :] * np.ones((n_dictionary, 1))
```
In [32]:
```def plot_scatter_MpVsComp(sparse_vector, my_sparse_code, title='MP', xlabel='MP', ylabel='COMP'):
fig = plt.figure(figsize=(16, 16))
a_min = np.min((sparse_vector.min(), my_sparse_code.min()))
a_max = np.max((sparse_vector.max(), my_sparse_code.max()))
ax.plot(np.array([a_min, a_max]), np.array([a_min, a_max]), 'k--', lw=2)
print(sparse_vector.shape, my_sparse_code.shape)
ax.scatter(sparse_vector.ravel(), my_sparse_code.ravel(), alpha=0.01)
ax.set_title(title)
ax.set_ylabel(ylabel)
ax.set_xlabel(xlabel)
#ax.set_xlim(0)
#ax.set_ylim(0)
ax.axis('equal')
return fig, ax

p_c = rescaling(corr, C=C, do_sym=do_sym)
fig, ax = plot_scatter_MpVsComp(p_c, quantile(P_cum_zeroeffect, p_c, stick),
title='transformation with flat pcum', xlabel='input', ylabel='quantile')
```
```(81792, 324) (81792, 324)
```

The relative mean squared error is low:

In [33]:
```print('Relative difference = ', np.sum( (p_c - quantile(P_cum_zeroeffect, p_c, stick))**2) / np.sum(p_c**2))
```
```Relative difference =  3.82686497427e-06
```

### COMP : using modulations¶

let's use this new `rescaling` function

In [34]:
```l0_sparseness = shl.l0_sparseness
def comp(data, dico, P_cum, C=C, do_sym=do_sym, verbose=0):
if verbose!=0: t0 = time.time()
n_samples, n_dictionary = data.shape[0], dico.shape[0]
sparse_code = np.zeros((n_samples, n_dictionary))
corr = (data @ dico.T)
Xcorr = (dico @ dico.T)
nb_quant = P_cum.shape[1]
stick = np.arange(n_dictionary)*nb_quant

for i_sample in range(n_samples):
c = corr[i_sample, :].copy()
if verbose!=0: ind_list=list()
for i_l0 in range(int(l0_sparseness)):
if P_cum is None:
q_i = rescaling(c, C=C, do_sym=do_sym)
else:
q_i = quantile(P_cum, rescaling(c, C=C, do_sym=do_sym), stick, do_fast=True)
ind  = np.argmax(q_i)
if verbose!=0: ind_list.append(ind)

c_ind = c[ind] / Xcorr[ind, ind]
sparse_code[i_sample, ind] += c_ind
c -= c_ind * Xcorr[ind, :]

if verbose!=0 and i_sample in range(2):
q_i = quantile(P_cum, rescaling(c, C=C, do_sym=do_sym), stick)
print(ind_list, [q_i[i] for i in ind_list], np.median(q_i), q_i.max(), [c[i] for i in ind_list], c.min(), c.max())
if verbose!=0:
duration = time.time()-t0
print('coding duration : {0}s'.format(duration))
return sparse_code

corr = (data_test @ dico_partial_learning.dictionary.T)
C = get_rescaling(corr, nb_quant=nb_quant, verbose=True)
sparse_code = comp(data_test, dico_partial_learning.dictionary, P_cum_zeroeffect, C=C, do_sym=do_sym, verbose=1)
```
```At indices (ranks)  [0, 51860, 103720, 155580, 207441, 259301, 311161, 363022, 414882, 466742, 518602, 570463, 622323, 674183, 726044, 777904, 829764, 881624, 933485, 985345, 1037205, 1089066, 1140926, 1192786, 1244646, 1296507, 1348367, 1400227, 1452088, 1503948, 1555808, 1607668, 1659529, 1711389, 1763249, 1815110, 1866970, 1918830, 1970690, 2022551, 2074411, 2126271, 2178132, 2229992, 2281852, 2333712, 2385573, 2437433, 2489293, 2541154, 2593014, 2644874, 2696734, 2748595, 2800455, 2852315, 2904176, 2956036, 3007896, 3059756, 3111617, 3163477, 3215337, 3267198, 3319058, 3370918, 3422778, 3474639, 3526499, 3578359, 3630220, 3682080, 3733940, 3785801, 3837661, 3889521, 3941381, 3993242, 4045102, 4096962, 4148823, 4200683, 4252543, 4304403, 4356264, 4408124, 4459984, 4511845, 4563705, 4615565, 4667425, 4719286, 4771146, 4823006, 4874867, 4926727, 4978587, 5030447, 5082308, 5134168, 5186028, 5237889, 5289749, 5341609, 5393469, 5445330, 5497190, 5549050, 5600911, 5652771, 5704631, 5756491, 5808352, 5860212, 5912072, 5963933, 6015793, 6067653, 6119513, 6171374, 6223234, 6275094, 6326955, 6378815, 6430675, 6482535, 6534396, 6586256, 6638116, 6689977, 6741837, 6793697, 6845557, 6897418, 6949278, 7001138, 7052999, 7104859, 7156719, 7208579, 7260440, 7312300, 7364160, 7416021, 7467881, 7519741, 7571602, 7623462, 7675322, 7727182, 7779043, 7830903, 7882763, 7934624, 7986484, 8038344, 8090204, 8142065, 8193925, 8245785, 8297646, 8349506, 8401366, 8453226, 8505087, 8556947, 8608807, 8660668, 8712528, 8764388, 8816248, 8868109, 8919969, 8971829, 9023690, 9075550, 9127410, 9179270, 9231131, 9282991, 9334851, 9386712, 9438572, 9490432, 9542292, 9594153, 9646013, 9697873, 9749734, 9801594, 9853454, 9905314, 9957175, 10009035, 10060895, 10112756, 10164616, 10216476, 10268336, 10320197, 10372057, 10423917, 10475778, 10527638, 10579498, 10631358, 10683219, 10735079, 10786939, 10838800, 10890660, 10942520, 10994380, 11046241, 11098101, 11149961, 11201822, 11253682, 11305542, 11357403, 11409263, 11461123, 11512983, 11564844, 11616704, 11668564, 11720425, 11772285, 11824145, 11876005, 11927866, 11979726, 12031586, 12083447, 12135307, 12187167, 12239027, 12290888, 12342748, 12394608, 12446469, 12498329, 12550189, 12602049, 12653910, 12705770, 12757630, 12809491, 12861351, 12913211, 12965071, 13016932, 13068792, 13120652, 13172513, 13224373, 13276233, 13328093, 13379954, 13431814, 13483674, 13535535, 13587395, 13639255, 13691115, 13742976, 13794836, 13846696, 13898557, 13950417, 14002277, 14054137, 14105998, 14157858, 14209718, 14261579, 14313439, 14365299, 14417159, 14469020, 14520880, 14572740, 14624601, 14676461, 14728321, 14780181, 14832042, 14883902, 14935762, 14987623, 15039483, 15091343, 15143204, 15195064, 15246924, 15298784, 15350645, 15402505, 15454365, 15506226, 15558086, 15609946, 15661806, 15713667, 15765527, 15817387, 15869248, 15921108, 15972968, 16024828, 16076689, 16128549, 16180409, 16232270, 16284130, 16335990, 16387850, 16439711, 16491571, 16543431, 16595292, 16647152, 16699012, 16750872, 16802733, 16854593, 16906453, 16958314, 17010174, 17062034, 17113894, 17165755, 17217615, 17269475, 17321336, 17373196, 17425056, 17476916, 17528777, 17580637, 17632497, 17684358, 17736218, 17788078, 17839938, 17891799, 17943659, 17995519, 18047380, 18099240, 18151100, 18202960, 18254821, 18306681, 18358541, 18410402, 18462262, 18514122, 18565982, 18617843, 18669703, 18721563, 18773424, 18825284, 18877144, 18929004, 18980865, 19032725, 19084585, 19136446, 19188306, 19240166, 19292027, 19343887, 19395747, 19447607, 19499468, 19551328, 19603188, 19655049, 19706909, 19758769, 19810629, 19862490, 19914350, 19966210, 20018071, 20069931, 20121791, 20173651, 20225512, 20277372, 20329232, 20381093, 20432953, 20484813, 20536673, 20588534, 20640394, 20692254, 20744115, 20795975, 20847835, 20899695, 20951556, 21003416, 21055276, 21107137, 21158997, 21210857, 21262717, 21314578, 21366438, 21418298, 21470159, 21522019, 21573879, 21625739, 21677600, 21729460, 21781320, 21833181, 21885041, 21936901, 21988761, 22040622, 22092482, 22144342, 22196203, 22248063, 22299923, 22351783, 22403644, 22455504, 22507364, 22559225, 22611085, 22662945, 22714806, 22766666, 22818526, 22870386, 22922247, 22974107, 23025967, 23077828, 23129688, 23181548, 23233408, 23285269, 23337129, 23388989, 23440850, 23492710, 23544570, 23596430, 23648291, 23700151, 23752011, 23803872, 23855732, 23907592, 23959452, 24011313, 24063173, 24115033, 24166894, 24218754, 24270614, 24322474, 24374335, 24426195, 24478055, 24529916, 24581776, 24633636, 24685496, 24737357, 24789217, 24841077, 24892938, 24944798, 24996658, 25048518, 25100379, 25152239, 25204099, 25255960, 25307820, 25359680, 25411540, 25463401, 25515261, 25567121, 25618982, 25670842, 25722702, 25774562, 25826423, 25878283, 25930143, 25982004, 26033864, 26085724, 26137584, 26189445, 26241305, 26293165, 26345026, 26396886, 26448746, 26500607]
the coefficients are  [-36.2 -6.93 -5.81 -5.19 -4.77 -4.44 -4.17 -3.95 -3.76 -3.59 -3.45 -3.31
-3.19 -3.08 -2.98 -2.89 -2.8 -2.72 -2.64 -2.57 -2.5 -2.43 -2.37 -2.31
-2.26 -2.2 -2.15 -2.1 -2.05 -2.01 -1.97 -1.92 -1.88 -1.84 -1.81 -1.77
-1.73 -1.7 -1.67 -1.63 -1.6 -1.57 -1.54 -1.51 -1.49 -1.46 -1.43 -1.41
-1.38 -1.36 -1.33 -1.31 -1.29 -1.26 -1.24 -1.22 -1.2 -1.18 -1.16 -1.14
-1.12 -1.1 -1.08 -1.07 -1.05 -1.03 -1.01 -0.997 -0.98 -0.964 -0.949 -0.933
-0.918 -0.903 -0.888 -0.874 -0.86 -0.846 -0.832 -0.819 -0.806 -0.793 -0.78
-0.768 -0.755 -0.743 -0.731 -0.72 -0.708 -0.697 -0.686 -0.675 -0.664
-0.653 -0.643 -0.632 -0.622 -0.612 -0.602 -0.592 -0.583 -0.573 -0.564
-0.555 -0.546 -0.537 -0.528 -0.519 -0.511 -0.502 -0.494 -0.486 -0.478
-0.47 -0.462 -0.454 -0.447 -0.439 -0.432 -0.424 -0.417 -0.41 -0.403 -0.396
-0.389 -0.382 -0.375 -0.369 -0.362 -0.356 -0.349 -0.343 -0.337 -0.331
-0.325 -0.319 -0.313 -0.307 -0.301 -0.295 -0.29 -0.284 -0.279 -0.273
-0.268 -0.263 -0.258 -0.252 -0.247 -0.242 -0.238 -0.233 -0.228 -0.223
-0.218 -0.214 -0.209 -0.205 -0.2 -0.196 -0.191 -0.187 -0.183 -0.179 -0.175
-0.17 -0.167 -0.163 -0.159 -0.155 -0.151 -0.147 -0.144 -0.14 -0.136 -0.133
-0.129 -0.126 -0.123 -0.119 -0.116 -0.113 -0.11 -0.107 -0.103 -0.1 -0.0975
-0.0946 -0.0918 -0.0891 -0.0864 -0.0837 -0.0811 -0.0785 -0.076 -0.0736
-0.0712 -0.0689 -0.0667 -0.0645 -0.0623 -0.0603 -0.0583 -0.0564 -0.0546
-0.0528 -0.0511 -0.0495 -0.0479 -0.0463 -0.0447 -0.0433 -0.0419 -0.0405
-0.0391 -0.0378 -0.0366 -0.0353 -0.0341 -0.033 -0.0318 -0.0307 -0.0296
-0.0286 -0.0276 -0.0265 -0.0256 -0.0246 -0.0236 -0.0226 -0.0217 -0.0208
-0.0199 -0.019 -0.0181 -0.0173 -0.0164 -0.0156 -0.0148 -0.014 -0.0132
-0.0124 -0.0116 -0.0108 -0.01 -0.00922 -0.00844 -0.00766 -0.00685 -0.00608
-0.00532 -0.00456 -0.00381 -0.00307 -0.00233 -0.00158 -0.000822 -7.28e-05
0.000709 0.00149 0.00227 0.00307 0.00386 0.00461 0.00537 0.00613 0.00693
0.00772 0.00851 0.00929 0.0101 0.0109 0.0117 0.0126 0.0134 0.0142 0.0151
0.016 0.0169 0.0178 0.0187 0.0196 0.0206 0.0215 0.0225 0.0235 0.0245
0.0255 0.0265 0.0275 0.0286 0.0297 0.0308 0.0321 0.0333 0.0345 0.0358
0.0371 0.0384 0.0398 0.0413 0.0428 0.0443 0.046 0.0477 0.0494 0.0511
0.0529 0.0547 0.0567 0.0587 0.0607 0.0628 0.065 0.0671 0.0693 0.0717
0.0741 0.0766 0.0792 0.0819 0.0846 0.0873 0.0901 0.093 0.0959 0.0989 0.102
0.105 0.108 0.111 0.115 0.118 0.121 0.125 0.128 0.132 0.135 0.139 0.143
0.146 0.15 0.154 0.158 0.162 0.166 0.17 0.174 0.178 0.183 0.187 0.191
0.195 0.2 0.204 0.209 0.214 0.218 0.223 0.228 0.233 0.238 0.243 0.248
0.253 0.258 0.263 0.269 0.274 0.28 0.285 0.291 0.296 0.302 0.308 0.314
0.32 0.326 0.332 0.338 0.344 0.351 0.357 0.363 0.37 0.377 0.384 0.39 0.397
0.405 0.412 0.419 0.426 0.434 0.441 0.449 0.457 0.464 0.472 0.48 0.489
0.497 0.505 0.514 0.522 0.531 0.54 0.549 0.558 0.567 0.577 0.586 0.596
0.606 0.615 0.626 0.636 0.646 0.657 0.667 0.678 0.689 0.701 0.712 0.723
0.735 0.747 0.759 0.772 0.784 0.797 0.81 0.823 0.836 0.85 0.864 0.878
0.893 0.907 0.922 0.937 0.952 0.968 0.984 1 1.02 1.03 1.05 1.07 1.09 1.11
1.12 1.14 1.16 1.18 1.2 1.22 1.25 1.27 1.29 1.31 1.34 1.36 1.39 1.41 1.44
1.46 1.49 1.52 1.55 1.58 1.61 1.64 1.67 1.7 1.74 1.77 1.81 1.85 1.89 1.93
1.97 2.01 2.06 2.1 2.15 2.2 2.26 2.31 2.37 2.43 2.5 2.57 2.64 2.71 2.8
2.88 2.98 3.08 3.19 3.31 3.44 3.59 3.75 3.94 4.16 4.42 4.74 5.17 5.8 6.91
36.6]
[112, 33, 278, 97, 101, 256, 30, 177, 36, 58, 59, 278, 88, 22, 278, 317] [0.66653387130965935, 0.65052856048959673, 0.59830775238067113, 0.62302444001996959, 0.64458959926104442, 0.64085691898403307, 0.61400716999153337, 0.57881440167231402, 0.0, 0.55064929939865515, 0.62791281020614631, 0.59830775238067113, 0.53729111115839245, 0.66534963023112492, 0.59830775238067113, 0.0] 0.0 0.726428524852 [0.14622114532348629, 0.11726767443482161, 0.04961471253073943, 0.075987077275044848, 0.10740013301351933, 0.10148169743045644, 0.065259558267194429, 0.034756512475519329, -0.026039702037990109, 0.019488148020597243, 0.082498085429622031, 0.04961471253073943, 0.013389142092147963, 0.14397244131880477, 0.04961471253073943, 0.0] -0.999707454888 0.287689288868
[200, 33, 73, 266, 6, 308, 278, 157, 137, 162, 194, 30, 290, 274, 212, 278] [0.52657777111988913, 0.0, 0.52231210787367488, 0.0, 0.53100666664042351, 0.50860769300083175, 0.0, 0.51572938871061569, 0.50472799517621458, 0.51739239907578116, 0.50876069676420965, 0.51956589106510176, 0.51830653729577003, 0.52029190910548584, 0.50996317355797349, 0.0] 0.0 0.550232245702 [0.0089410744039219933, -0.013120664016748315, 0.0072309348418384935, -0.00013820264737959961, 0.010729581615201201, 0.0017983779666867076, 0.0, 0.0046274909864738086, 0.00024940207037818778, 0.0052722273582512689, 0.001859008757976976, 0.006117370455199618, 0.0056274915879661741, 0.0064112040165253299, 0.0023374950692405774, 0.0] -0.07863481717 0.0192883783941
coding duration : 324.98530888557434s
```
In [35]:
```%%timeit
corr = (data_test @ dico_partial_learning.dictionary.T)
C = get_rescaling(corr, nb_quant=nb_quant, verbose=False)
sparse_code = comp(data_test, dico_partial_learning.dictionary, P_cum_zeroeffect, C=C, do_sym=do_sym, verbose=0)
```
```3min 24s ± 29.3 s per loop (mean ± std. dev. of 7 runs, 1 loop each)
```
In [36]:
```%%timeit
corr = (data_test @ dico_partial_learning.dictionary.T)
C = get_rescaling(corr, nb_quant=nb_quant, verbose=False)
sparse_code = comp(data_test, dico_partial_learning.dictionary, P_cum_zeroeffect, C=C, do_sym=do_sym, verbose=0)
```
```3min 1s ± 6.55 s per loop (mean ± std. dev. of 7 runs, 1 loop each)
```

### testing that COMP with fixed Pcum is equivalent to MP¶

In [37]:
```print(dico_partial_learning.P_cum)
```
```None
```
In [38]:
```from shl_scripts.shl_learn import get_P_cum
from shl_scripts.shl_encode import get_rescaling, rescaling
corr = (data_test @ dico_partial_learning.dictionary.T)
C = get_rescaling(corr, nb_quant=nb_quant, verbose=True)
p_c = rescaling(corr, C=C, do_sym=do_sym)
```
In [39]:
```from shl_scripts.shl_encode import get_rescaling, mp
l0_sparseness = 50

sparse_code_mp = mp(data_test, dico_partial_learning.dictionary, P_cum=None, C=5., do_sym=do_sym, l0_sparseness=l0_sparseness, verbose=1)
sparse_code_comp = mp(data_test, dico_partial_learning.dictionary, P_cum=P_cum_zeroeffect, C=5., do_sym=do_sym, l0_sparseness=l0_sparseness, verbose=1)

#print('Difference = ', np.sum((sparse_code_mp - sparse_code_comp)**2))
#print('Relative difference = ', np.mean((sparse_code_mp - sparse_code_comp)**2))
#print('Variance = ', np.mean((sparse_code_mp)**2))
print('Relative difference = ', np.sum((sparse_code_mp - sparse_code_comp)**2)/np.sum((sparse_code_mp)**2))

fig, ax = plot_proba_histogram(sparse_code_mp)
fig.show()
fig, ax = plot_proba_histogram(sparse_code_comp)
```
```coding duration : 90.52337574958801
coding duration : 346.8867189884186
Relative difference =  0.0
```
In [40]:
```n_samples, nb_filter = sparse_code_mp.shape

from shl_scripts.shl_encode import get_rescaling, mp

corr = (data_test @ dico_partial_learning.dictionary.T)
C_vec = get_rescaling(corr, nb_quant=nb_quant, verbose=True)

P_cum = np.vstack((P_cum_zeroeffect, C_vec))

print('Classical MP')
sparse_code_mp = mp(data_test, dico_partial_learning.dictionary, P_cum=None, C=np.inf, do_sym=do_sym, l0_sparseness=l0_sparseness, verbose=1)
print('Homeostatic MP')
sparse_code_comp = mp(data_test, dico_partial_learning.dictionary, P_cum=P_cum, C=0., do_sym=do_sym, l0_sparseness=l0_sparseness, verbose=1)

print('Relative difference = ', np.sum((sparse_code_mp - sparse_code_comp)**2)/np.sum((sparse_code_mp)**2))
```
```Classical MP
coding duration : 89.88441371917725
Homeostatic MP
coding duration : 536.8180718421936
Relative difference =  0.122805037453
```
In [41]:
```fig, ax = plot_proba_histogram(sparse_code_mp)
fig.show()
fig, ax = plot_proba_histogram(sparse_code_comp)
```
In [42]:
```fig, ax = plot_scatter_MpVsComp(sparse_code_mp, sparse_code_comp)
```
```(81792, 324) (81792, 324)
```