<jupyter_start><jupyter_text>x and y data<jupyter_code>error_probability_list=np.array(error_probability_list)

pc=0.280
nu=1.3
rescaled_x=[(error_probability_list-pc)*9**(1/nu),(error_probability_list-pc)*13**(1/nu),(error_probability_list-pc)*17**(1/nu),(error_probability_list-pc)*21**(1/nu)]
rescaled_x=np.reshape(rescaled_x,[4,11])

# fit a second degree polynomial to the economic data
from numpy import arange
from pandas import read_csv
from scipy.optimize import curve_fit
from matplotlib import pyplot
 
# Pfail = A + Bx + Cx2 + DL^(−1/µ), where A, B and C are expansion coefficients, D is the coefficient of the 
# non-universal term, and x = (p − pc0)L^(1/ν0)

# define the true objective function
def objective(x, a, b, c):
	return a * x + b * x**2 + c
 
# load the dataset
# choose the input and output variables
x, y = rescaled_x[:,4:9], Pfail[:,4:9]
# curve fit
popt0, _ = curve_fit(objective, x[0], y[0])
popt1, _ = curve_fit(objective, x[1], y[1])
popt2, _ = curve_fit(objective, x[2], y[2])
popt3, _ = curve_fit(objective, x[3], y[3])
# summarize the parameter values
a, b, c = popt0
print('y = %.5f * x + %.5f * x^2 + %.5f' % (a, b, c))

# plot input vs output
pyplot.scatter(x[0], y[1])

# define a sequence of inputs between the smallest and largest known inputs
x_line = x[0]

# calculate the output for the range
y_line = objective(x_line, a, b, c)
# create a line plot for the mapping function
pyplot.plot(x_line, y_line, '--', color='red')
pyplot.show()<jupyter_output>y = 0.53432 * x + -0.17348 * x^2 + 0.27377<jupyter_text>Using the pyfssa package<jupyter_code>import itertools
import fssa
from cycler import cycler
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns

# Define colors
palette = sns.cubehelix_palette(
    n_colors=10, start=2.0, rot=0.35, gamma=1.0, hue=1.0, light=0.6, dark=0.2,
)
sns.palplot(palette)

#we use the colors below instead of the above
palette = plt.cm.Spectral(np.linspace(0,10,50))

#choose various values of the critical exponents 
pX_c = np.tile(0.28, (3, 3))
nu = np.tile(1.45, (3, 3))
zeta = np.tile(0, (3, 3))
pX_c[0, :] = [0.28,0.28, 0.28]
nu[1, :] = [1.05, 1.1, 1.15]
zeta[2, :] = [0, 0,0]
std_list_avg=np.array(Pfail_err)
pL_list_avg=Pfail
pX_list=error_probability_list

# re-scale data (manually)
scaled_data = list()
quality = list()
for i in range(3):
    my_scaled_data = list()
    my_quality = list()
    for j in range(3):
        my_scaled_data.append(
            fssa.scaledata(
                L_list,  pX_list, pL_list_avg, std_list_avg,
                pX_c[i, j], nu[i, j], zeta[i, j]
            )
        )
        my_quality.append(fssa.quality(*my_scaled_data[-1]))
    scaled_data.append(my_scaled_data)
    quality.append(my_quality)
    
        
# plot manually re-scaled data
fig, axes = plt.subplots(
    nrows=3, ncols=3, squeeze=True,
    #figsize=(8, 7),
    sharex=True, sharey=True,
)

for (i, j) in itertools.product(range(3), range(3)):
    ax = axes[i, j]
    ax.set_prop_cycle(cycler('color', palette))
    my_scaled_data = scaled_data[i][j]
    for l_index, l in enumerate(L_list):
        ax.plot(
            my_scaled_data.x[l_index, :], my_scaled_data.y[l_index, :],
            '.',
            label=r'${}$'.format(l),
            rasterized=True,
        )
    ax.set_xbound(-0.4,0.7)
    if i == 0:
        ax.set_title(
            r'$pX_c = {}$'.format(pX_c[i, j]),
            position=(0.15, 0.65),
        )
    elif i == 1:
        ax.set_title(
            r'$\nu = {}$'.format(nu[i, j]),
            position=(0.25, 0.65),
        )
    elif i == 2:
        ax.set_title(
            r'$\zeta = {}$'.format(zeta[i, j]),
            position=(0.25, 0.65),
        )
    if i == 2:
        ax.set_xlabel(r'$x$')
        ax.set_xticks([0,0.2,0.3,0.4,0.5,0.6])
    if j == 0:
        ax.set_yticks([0, 1, 2, 3, 4, 5])
    ax.text(
        0.1, 0.5,
        r'$S={:.1f}$'.format(quality[i][j]),
        transform=ax.transAxes,
    )

plt.show()
ret = fssa.autoscale(L_list, pX_list, pL_list_avg, std_list_avg, 0.5, 1.4, 0)
auto_scaled_data = fssa.scaledata(L_list, pX_list, pL_list_avg, std_list_avg,  ret.rho, ret.nu, ret.zeta)

# critical exponents and errors, quality of data collapse
print(ret.rho, ret.drho)
print(ret.nu, ret.dnu)
print(ret.zeta, ret.dzeta)
print(ret.fun)

fig, ax = plt.subplots()
ax.set_prop_cycle(cycler('color', palette))
ax.plot(
    auto_scaled_data.x.T, auto_scaled_data.y.T,
    '.',
)
ax.set_xbound(-4, 2)
ax.set_xlabel(r'$x$')
plt.show()<jupyter_output>/opt/sagemath-9.2/local/lib/python3.7/site-packages/fssa/fssa.py:165: RuntimeWarning: rho_c is out of range
  warnings.warn("rho_c is out of range", RuntimeWarning)
/opt/sagemath-9.2/local/lib/python3.7/site-packages/fssa/fssa.py:505: RuntimeWarning: invalid value encountered in sqrt
  return np.sqrt(np.diag(varco)), varco
{'fraction_non_alphanumeric': 0.16036960985626283, 'fraction_numerical': 0.04188911704312115, 'mean_word_length': 3.3106194690265487, 'pattern_counts': {'":': 0, '<': 7, '<?xml version=': 0, '>': 7, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '644474', 'license': 'BSD-3-Clause', 'n_tokens_mistral': 2111, 'n_tokens_neox': 2017, 'n_words': 436, 'path': 'experiments/find-threshold-randXZYZ/scaling_bias10.ipynb', 'repo_name': 'dua-arpit/qecsim', 'size': '4870'}