<jupyter_start><jupyter_text>Get data<jupyter_code>random.seed(seed)
data_manager = data_loader.AMPDataManager(
    '../data/unlabelled_positive.csv',
    '../data/unlabelled_negative.csv',
    min_len=MIN_LENGTH,
    max_len=MAX_LENGTH)

amp_x, amp_y = data_manager.get_merged_data()
amp_x_train, amp_x_test, amp_y_train, amp_y_test = train_test_split(amp_x, amp_y, test_size=0.1, random_state=36)
amp_x_train, amp_x_val, amp_y_train, amp_y_val = train_test_split(amp_x_train, amp_y_train, test_size=0.2, random_state=36)
# Restrict the length
ecoli_df = pd.read_csv('../data/mic_data.csv')
mask = (ecoli_df['sequence'].str.len() <= MAX_LENGTH) & (ecoli_df['sequence'].str.len() >= MIN_LENGTH)
ecoli_df = ecoli_df.loc[mask]
mic_x = pad(to_one_hot(ecoli_df['sequence']))
mic_y = ecoli_df.value

mic_x_train, mic_x_test, mic_y_train, mic_y_test = train_test_split(mic_x, mic_y, test_size=0.1, random_state=36)
mic_x_train, mic_x_val, mic_y_train, mic_y_val = train_test_split(mic_x_train, mic_y_train, test_size=0.2, random_state=36)
n=64

pos = np.vstack([amp_x_test[amp_y_test == 1], mic_x_test[mic_y_test < 1.5]])
neg = np.vstack([amp_x_test[amp_y_test == 0], mic_x_test[mic_y_test > 1.5]])

pos_amp = amp_classifier_model.predict(pos, verbose=1).reshape(len(pos))
neg_mic = mic_classifier_model.predict(neg, verbose=1).reshape(len(neg))
neg_amp = amp_classifier_model.predict(neg, verbose=1).reshape(len(neg))
pos_mic = mic_classifier_model.predict(pos, verbose=1).reshape(len(pos))

pos = np.vstack([pos] * n).reshape(-1, 25)
pos_amp = np.vstack([pos_amp] * n).reshape(-1, 1)
pos_mic = np.vstack([pos_mic] * n).reshape(-1, 1)
neg = np.vstack([neg] * n).reshape(-1, 25)
neg_amp = np.vstack([neg_amp] * n).reshape(-1, 1)
neg_mic = np.vstack([neg_mic] * n).reshape(-1, 1)
def get_z_sigma(x):
    inputs = layers.Input(shape=(25,))
    z_mean, z_sigma, z = AMPMaster.encoder.output_tensor(inputs)
    temp_encoder = m.Model(inputs, [z_mean, z_sigma, z])
    z_mean, z_sigma, z = temp_encoder.predict(x) 
    return np.exp(z_sigma / 2)
def improve(seq, amp, mic, z_sigma, temp=0.0):
    z = encoder_model.predict(seq, batch_size=5000)
    noise = np.random.normal(loc=0, scale=temp*z_sigma, size=z.shape)
    encoded = z + noise
    conditioned = np.hstack([
        encoded,
        np.ones((len(seq), 1)),
        np.ones((len(seq), 1)),
    ])
    decoded = decoder_model.predict(conditioned, batch_size=5000)
    new_peptides = np.argmax(decoded, axis=2)
    new_amp = amp_classifier_model.predict(new_peptides, batch_size=5000)
    new_mic = mic_classifier_model.predict(new_peptides, batch_size=5000)                                                                             
    
    # RELATIVE
    rel_better = new_amp > amp.reshape(-1, 1)
    rel_better = rel_better & (new_mic > mic.reshape(-1, 1))
    rel_better = np.logical_or.reduce(rel_better, axis=1)
    rel_improved = new_peptides[np.where(rel_better), :].reshape(-1, 25)
    before_rel_improve = seq[np.where(rel_better), :].reshape(-1, 25)
    
    # ABSOLUTE
    abs_better = new_amp >= 0.8
    abs_better = abs_better & (new_mic > 0.5)
    abs_better = np.logical_or.reduce(abs_better, axis=1)
    abs_improved = new_peptides[np.where(abs_better), :].reshape(-1, 25)
    before_abs_improve = seq[np.where(abs_better), :].reshape(-1, 25)
    
    return {
        'new_peptides': new_peptides,
        'rel_improved': rel_improved,
        'abs_improved': abs_improved,
        'before_rel_improve': before_rel_improve,
        'before_abs_improve': before_abs_improve, 
        'new_amp': new_amp,
        'new_mic': new_mic,
        }                                                                  
tau = backend.variable(1.0, name="temperature")
model = 'HydrAMP'
AMPMaster = bms.load_model(f'../models/{model}/{best_epochs[model]}')
encoder_model =  AMPMaster.encoder(input_to_encoder)
AMPMaster.decoder.activation.temperature = tau
decoder_model = AMPMaster.decoder(input_to_decoder)<jupyter_output>tracking <tf.Variable 'temperature:0' shape=() dtype=float32, numpy=0.1> temperature
tracking <tf.Variable 'temperature:0' shape=() dtype=float32, numpy=1.0> temperature<jupyter_text>Positives, negatives<jupyter_code>temps =  np.linspace(0,5,num=21)
pos_analogues = []
neg_analogues = []
pos_z_sigma = get_z_sigma(pos)
neg_z_sigma = get_z_sigma(neg)
for temp in tqdm(temps):
    pos_analogues.append(np.unique(improve(pos, pos_amp, pos_mic, pos_z_sigma, temp=temp)['rel_improved'], axis=0))
    neg_analogues.append(np.unique(improve(neg, neg_amp, neg_mic, neg_z_sigma, temp=temp)['abs_improved'], axis=0))
no_pos_analogues = [len(x) for x in pos_analogues]
no_neg_analogues = [len(x) for x in neg_analogues]
fig, (ax1, ax2) = plt.subplots(
    nrows=1,
    ncols=2,
    figsize=(12, 4),
    dpi=90,
)

ax1.bar(temps, no_pos_analogues, width=0.25)
ax1.set_title('Positives: relative improvement')
ax2.bar(temps, no_neg_analogues, width=0.25)
ax2.set_title('Negatives: absolute improvement')

plt.show()
neg
neg0 = improve(neg, neg_amp, neg_mic, neg_z_sigma, temp=0)
neg2 = improve(neg, neg_amp, neg_mic, neg_z_sigma, temp=2)
neg5 = improve(neg, neg_amp, neg_mic, neg_z_sigma, temp=5)

dump(neg0, f'../results/z_sigma_improvement_0_NegTest.joblib')
dump(neg2, f'../results/z_sigma_improvement_2_NegTest.joblib')
dump(neg5, f'../results/z_sigma_improvement_5_NegTest.joblib')
len(neg_amp)
pos0 = improve(pos, pos_amp, pos_mic, pos_z_sigma, temp=0)
pos2 = improve(pos, pos_amp, pos_mic, pos_z_sigma, temp=2)
pos5 = improve(pos, pos_amp, pos_mic, pos_z_sigma, temp=5)

dump(pos0, f'../results/z_sigma_improvement_0_PosTest.joblib')
dump(pos2, f'../results/z_sigma_improvement_2_PosTest.joblib')
dump(pos5, f'../results/z_sigma_improvement_5_PosTest.joblib')<jupyter_output><empty_output><jupyter_text>Pexiganan, tenporin, CAMEL<jupyter_code>def levenshteinDistanceDP(token1, token2):
    distances = np.zeros((len(token1) + 1, len(token2) + 1))

    for t1 in range(len(token1) + 1):
        distances[t1][0] = t1

    for t2 in range(len(token2) + 1):
        distances[0][t2] = t2
        
    a = 0
    b = 0
    c = 0
    
    for t1 in range(1, len(token1) + 1):
        for t2 in range(1, len(token2) + 1):
            if (token1[t1-1] == token2[t2-1]):
                distances[t1][t2] = distances[t1 - 1][t2 - 1]
            else:
                a = distances[t1][t2 - 1]
                b = distances[t1 - 1][t2]
                c = distances[t1 - 1][t2 - 1]
                
                if (a <= b and a <= c):
                    distances[t1][t2] = a + 1
                elif (b <= a and b <= c):
                    distances[t1][t2] = b + 1
                else:
                    distances[t1][t2] = c + 1

    return distances[len(token1)][len(token2)]
def get_unique(x):
    if x.shape[0] > 1:
        return np.unique(x, axis=0)
    else:
        return x
peptides = [
    'GIGKFLKKAKKFGKAFVKILKK', #pex
    'FLPLIGRVLSGIL', #temporin
    'KWKLFKKIGAVLKVL'  #camel
]

n = 10000
temps =  np.linspace(0,5,num=21)
no_analogues = []
levenshtein = []

for peptide in tqdm(peptides):
    pep = pad(to_one_hot([peptide]))
    amp = amp_classifier_model.predict(pep)
    mic = mic_classifier_model.predict(pep)
    z_sigma = get_z_sigma(pep)
    
    pep = np.vstack([pep] * n).reshape(-1, 25)
    amp = np.vstack([amp] * n)
    mic = np.vstack([mic] * n)
    
    analogues = []
    ld = []
    
    for temp in temps:
        analogue_batch = get_unique(improve(pep, amp, mic, z_sigma, temp=temp)['abs_improved'])
        analogues.append(analogue_batch)
        ld.append([levenshteinDistanceDP(peptide, translate_peptide(x)) for x in analogue_batch])

    no_analogues.append([len(x) for x in analogues])
    levenshtein.append(ld)

fig, (ax1, ax2) = plt.subplots(
    nrows=1,
    ncols=2,
    figsize=(12, 4),
    dpi=90,
)

# ax1.bar(temps, no_analogues[0], width=0.25)
ax1.bar(temps, np.log(no_analogues[0]), width=0.25)

ax1.set_title('Pexiganan: absolute improvement')

flierprops = dict(marker='o', markersize=1,
                  linestyle='none')
ax2.boxplot(levenshtein[0], sym='+', flierprops=flierprops)
ax2.set_xticklabels(temps, rotation=90)
ax2.set_title('Levenshtein distance')

plt.show()
fig, (ax1, ax2) = plt.subplots(
    nrows=1,
    ncols=2,
    figsize=(12, 4),
    dpi=90,
)

# ax1.bar(temps, no_analogues[1], width=0.25)
ax1.set_xlim(-0.25,10.5)
ax1.bar(temps, np.log(no_analogues[1]), width=0.25)
ax1.set_title('Temporin: absolute improvement')

flierprops = dict(marker='o', markersize=1,
                  linestyle='none')
ax2.boxplot(levenshtein[1], sym='+', flierprops=flierprops)
ax2.set_xticklabels(temps, rotation=90)
ax2.set_title('Levenshtein distance')
ax2.set_ylim(0, 25)

plt.show()
fig, (ax1, ax2) = plt.subplots(
    nrows=1,
    ncols=2,
    figsize=(12, 4),
    dpi=90,
)

# ax1.bar(temps, no_analogues[2], width=0.25)
ax1.bar(temps, np.log(no_analogues[2]), width=0.25)

ax1.set_title('Camel: absolute improvement')

flierprops = dict(marker='o', markersize=1,
                  linestyle='none')
ax2.boxplot(levenshtein[2], sym='+', flierprops=flierprops)
ax2.set_xticklabels(temps, rotation=90)
ax2.set_title('Levenshtein distance')

plt.show()
temps = [0, 2, 5]
no_analogues = []
levenshtein = []

for peptide in tqdm(peptides):
    pep = pad(to_one_hot([peptide]))
    amp = amp_classifier_model.predict(pep)
    mic = mic_classifier_model.predict(pep)
    z_sigma = get_z_sigma(pep)
    
    pep = np.vstack([pep] * n).reshape(-1, 25)
    amp = np.vstack([amp] * n)
    mic = np.vstack([mic] * n)
    
    analogues = []
    ld = []
    
    for temp in temps:
        analogue_batch = get_unique(improve(pep, amp, mic, z_sigma, temp=temp)['abs_improved'])
        analogues.append(analogue_batch)
        ld.append([levenshteinDistanceDP(peptide, translate_peptide(x)) for x in analogue_batch])

    no_analogues.append([len(x) for x in analogues])
    levenshtein.append(ld)

pex_temp_camel = {'no_analogues': no_analogues, "ld": levenshtein}
dump(pex_temp_camel, f'../results/PexTempCamel.joblib')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1486848635235732, 'fraction_numerical': 0.03225806451612903, 'mean_word_length': 3.7866983372921617, 'pattern_counts': {'":': 1, '<': 20, '<?xml version=': 0, '>': 21, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 6, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '301802', 'license': 'MIT', 'n_tokens_mistral': 4306, 'n_tokens_neox': 4027, 'n_words': 657, 'path': 'scripts/z_sigma_improvement.ipynb', 'repo_name': 'szczurek-lab/hydramp', 'size': '10075'}