analysis.org

Performance figures

from utils import paired_positions
import subprocess
import matplotlib.pyplot as plt
from matplotlib import rcParams, artist
from matplotlib.patches import Patch
from numpy import mean, std, median, linspace, arange
from scipy.stats import kde
from RNA import fold_compound, bp_distance
from utility.utils_analysis import get_loop_content, read_csv, read_true_struct
plt.rcParams["font.family"] = "serif"
plt.rcParams["font.size"] = 12
plt.rcParams["figure.figsize"] = 7, 4

# fft_pred = read_csv("./benchmark_results/fft_100n_50ms_best_nrj_scores.csv")
fft_pred = read_csv("./benchmark_results/fft_200n_200ms_scores.csv")
# fftb_pred = read_csv("./benchmark_results/fft_2_scores.csv")
fftb_pred = read_csv("./benchmark_results/fft_100n_50ms_scores.csv")
mx_pred = read_csv("./benchmark_results/mxfold_scores.csv")
vrna_pred = read_csv("./benchmark_results/mfe_scores.csv")
true_str = read_true_struct()
fft_all_v, mfe_all_v, fftb_all_v, mx_all_v = [], [], [], []
fft_to_plot, mfe_to_plot, fftb_to_plot, mx_to_plot = {}, {}, {}, {}

nb = 0
for seq, (struct, name) in true_str.items():
    seq_comp = fold_compound(seq)
    true_loop_content = get_loop_content(struct)

    true_nb_bp = struct.count("(")
    true_nrj = seq_comp.eval_structure(struct)
    len_seq, mfe_st, mfe_nrj, mfe_nbp, mfe_pvv, mfe_sens = vrna_pred[seq]

    # if seq in fft_pred and seq in fftb_pred and seq in mx_pred and len(seq) > 100 and len(seq) < 150:
    # if seq in fft_pred and seq in fftb_pred and seq in mx_pred and fft_pred[seq][1] != "."*len(seq):
    if seq in fft_pred and seq in fftb_pred and seq in mx_pred:
        len_seq, fft_st, fft_nrj, fft_nbp, fft_pvv, fft_sens = fft_pred[seq]
        len_seq, fftb_st, fftb_nrj, fftb_nbp, fftb_pvv, fftb_sens = fftb_pred[seq]
        len_seq, mx_st, mx_nrj, mx_nbp, mx_pvv, mx_sens = mx_pred[seq]

        save_fft, save_mfe, save_fftb, save_mx = fft_sens, mfe_sens, fftb_sens, mx_sens
        fft_all_v += [save_fft]
        mfe_all_v += [save_mfe]
        fftb_all_v += [save_fftb]
        mx_all_v += [save_mx]

        # if mfe_pvv <= 30.0:
        #     print(name)
        #     print(seq)
        #     print(fft_pvv, mfe_pvv, mx_pvv)
        #     print(fft_st)
        #     # print(mfe_st)
        #     # print(mx_st)
        #     print(struct)
        #     print()

        if len_seq in fft_to_plot:
            fft_to_plot[len_seq] += [save_fft]
            mfe_to_plot[len_seq] += [save_mfe]
            fftb_to_plot[len_seq] += [save_fftb]
            mx_to_plot[len_seq] += [save_mx]
        else:
            fft_to_plot[len_seq] = [save_fft]
            mfe_to_plot[len_seq] = [save_mfe]
            fftb_to_plot[len_seq] = [save_fftb]
            mx_to_plot[len_seq] = [save_mx]


list_len = list(fft_to_plot.keys())
list_len.sort()
nb_el = len(list_len)

fft_mean = mean([mean(fft_to_plot[len_]) for len_ in list_len if len_ > 0])
mfe_mean = mean([mean(mfe_to_plot[len_]) for len_ in list_len if len_ > 0])
fftb_mean = mean([mean(fftb_to_plot[len_]) for len_ in list_len if len_ > 0])
mx_mean = mean([mean(mx_to_plot[len_]) for len_ in list_len if len_ > 0])

slide_mean_fft = [mean([mean(fft_to_plot[l]) for l in list_len[max(0, i-10):min(nb_el, i+10)]]) for i, len_ in enumerate(list_len)]
slide_mean_mfe = [mean([mean(mfe_to_plot[l]) for l in list_len[max(0, i-10):min(nb_el, i+10)]]) for i, len_ in enumerate(list_len)]
slide_mean_fftb = [mean([mean(fftb_to_plot[l]) for l in list_len[max(0, i-10):min(nb_el, i+10)]]) for i, len_ in enumerate(list_len)]
slide_mean_mx = [mean([mean(mx_to_plot[l]) for l in list_len[max(0, i-10):min(nb_el, i+10)]]) for i, len_ in enumerate(list_len)]

bpros = dict(color="orange", linewidth=0.2)
fig = plt.figure(1)

left, width = 0.10, 0.69
bottom, height = 0.18, 0.79
spacing = 0.000
rect_scatter = [left, bottom, width, height]
rect_histy = [left + width + spacing, bottom, 0.2, height]
fft_f = fig.add_axes(rect_scatter)

ax_histy = fig.add_axes(rect_histy, sharey=fft_f)
ax_histy.set_ylim([0, 100.0])
xs_scr = linspace(min(fft_all_v), max(fft_all_v), 50)
dens_fft_scr = kde.gaussian_kde(fft_all_v)
dens_mfe_scr = kde.gaussian_kde(mfe_all_v)
dens_fftb_scr = kde.gaussian_kde(fftb_all_v)
dens_mx_scr = kde.gaussian_kde(mx_all_v)
ax_histy.plot(dens_fft_scr(xs_scr), xs_scr, c="deepskyblue", linewidth=0.8)
ax_histy.plot(dens_mfe_scr(xs_scr), xs_scr, c="orangered", linewidth=0.8)
ax_histy.plot(dens_fftb_scr(xs_scr), xs_scr, c="green", linewidth=0.8)
ax_histy.plot(dens_mx_scr(xs_scr), xs_scr, c="blueviolet", linewidth=0.8)
ax_histy.tick_params(axis="y", labelleft=False, size=0)
ax_histy.tick_params(axis="x", labelbottom=False, size=0)
ax_histy.grid(True, color="grey",linestyle="--", linewidth=0.2)

fft_f.set_xlim([0, nb_el]); fft_f.set_ylim([0, 100.0])
fft_f.grid(True, color="grey",linestyle="--", linewidth=0.2)

fft_f.errorbar([i+0.0 for i in range(nb_el)], [mean(fft_to_plot[len_]) for len_ in list_len],
               yerr=[std(fft_to_plot[len_])/4 for len_ in list_len], elinewidth=0.8, fmt=".",
               linewidth=0.2, ms=3, c="skyblue", ecolor="deepskyblue", alpha=0.2)
fft_f.errorbar([i+0.0 for i in range(nb_el)], [mean(mfe_to_plot[len_]) for len_ in list_len],
               yerr=[std(mfe_to_plot[len_])/4 for len_ in list_len], elinewidth=0.8, fmt=".",
               linewidth=0.2, ms=3, c="coral", ecolor="orangered", alpha=0.2)
fft_f.errorbar([i+0.0 for i in range(nb_el)], [mean(fftb_to_plot[len_]) for len_ in list_len],
               yerr=[std(fftb_to_plot[len_])/4 for len_ in list_len], elinewidth=0.8, fmt=".",
               linewidth=0.2, ms=3, c="palegreen", ecolor="palegreen", alpha=0.2)
fft_f.errorbar([i+0.0 for i in range(nb_el)], [mean(mx_to_plot[len_]) for len_ in list_len],
               yerr=[std(mx_to_plot[len_])/4 for len_ in list_len], elinewidth=0.8, fmt=".",
               linewidth=0.2, ms=3, c="blueviolet", ecolor="blueviolet", alpha=0.2)

fft_f.scatter([i+0.0 for i in range(nb_el)], [mean(fft_to_plot[len_]) for len_ in list_len], c="deepskyblue", s=0.5)
fft_f.scatter([i+0.0 for i in range(nb_el)], [mean(mfe_to_plot[len_]) for len_ in list_len], c="orangered", s=0.5)
fft_f.scatter([i+0.0 for i in range(nb_el)], [mean(fftb_to_plot[len_]) for len_ in list_len], c="green", s=0.5)
fft_f.scatter([i+0.0 for i in range(nb_el)], [mean(mx_to_plot[len_]) for len_ in list_len], c="blueviolet", s=0.5)

fft_f.plot(range(nb_el), slide_mean_fft, linestyle="--", c="deepskyblue", alpha = 0.8)
fft_f.plot(range(nb_el), slide_mean_mfe, linestyle="--", c="orangered", alpha = 0.8)
fft_f.plot(range(nb_el), slide_mean_fftb, linestyle="--", c="green", alpha = 0.8)
fft_f.plot(range(nb_el), slide_mean_mx, linestyle="--", c="blueviolet", alpha = 0.8)
# fft_f.set_xlabel("Length", fontsize=15)
fft_f.tick_params(axis="x", labelbottom=False)
fft_f.set_ylabel("PPV", fontsize=15)

fft_l = Patch(color="skyblue", label="FFT")
mfe_l = Patch(color="tomato", label="VRNA")
fftb_l = Patch(color="green", label="Line")
mx_l = Patch(color="blueviolet", label="ML")

fft_f.xaxis.set_ticklabels(list_len[::20], rotation=45)
fft_f.xaxis.set_ticks([i for i in range(nb_el)][::20])

name_fft = "RAFFT"
name_mfe = "MFE"
name_fftb = "RAFFT*"
name_mx = "ML"

fft_f.legend((fft_l,fftb_l,mfe_l,mx_l), (f"{name_fft}", f"{name_fftb}",
                                    f"{name_mfe}", f"{name_mx}"), loc="lower center", ncol=4, fontsize=10)
# plt.savefig("img/fold_perf_ppv_50.png", dpi=300, transparent=True)
print(f"{fft_mean:.1f} {fftb_mean:.1f} {mfe_mean:.1f} {mx_mean:.1f}")
plt.show()

def function(elems):
    """something
    """
    return 0

print(function(0))

Resutls statistical significance

from utils import paired_positions
from numpy import mean, std, median, linspace, arange
from scipy.stats import kde, ttest_ind
from RNA import fold_compound, bp_distance
from utility.utils_analysis import get_loop_content, read_csv, read_true_struct
from math import isnan

fft_pred = read_csv("./benchmark_results/fft_100n_50ms_best_nrj_scores.csv")
fftb_pred = read_csv("./benchmark_results/fft_100n_50ms_scores.csv")
mx_pred = read_csv("./benchmark_results/mxfold_scores.csv")
vrna_pred = read_csv("./benchmark_results/mfe_scores.csv")
true_str = read_true_struct()
fft_all_v, mfe_all_v, fftb_all_v, mx_all_v = [], [], [], []
fft_to_plot, mfe_to_plot, fftb_to_plot, mx_to_plot = {}, {}, {}, {}

nb = 0
for seq, (struct, name) in true_str.items():
    seq_comp = fold_compound(seq)
    true_loop_content = get_loop_content(struct)

    true_nb_bp = struct.count("(")
    true_nrj = seq_comp.eval_structure(struct)

    if seq in fft_pred and seq in fftb_pred and seq in mx_pred:
        len_seq, mfe_st, mfe_nrj, mfe_nbp, mfe_pvv, mfe_sens = vrna_pred[seq]
        len_seq, fft_st, fft_nrj, fft_nbp, fft_pvv, fft_sens = fft_pred[seq]
        len_seq, fftb_st, fftb_nrj, fftb_nbp, fftb_pvv, fftb_sens = fftb_pred[seq]
        len_seq, mx_st, mx_nrj, mx_nbp, mx_pvv, mx_sens = mx_pred[seq]

        save_fft, save_mfe, save_fftb, save_mx = fft_pvv, mfe_pvv, fftb_pvv, mx_pvv
        fft_all_v += [save_fft]
        mfe_all_v += [save_mfe]
        fftb_all_v += [save_fftb]
        mx_all_v += [save_mx]

        if len_seq in fft_to_plot:
            fft_to_plot[len_seq] += [save_fft]
            mfe_to_plot[len_seq] += [save_mfe]
            fftb_to_plot[len_seq] += [save_fftb]
            mx_to_plot[len_seq] += [save_mx]
        else:
            fft_to_plot[len_seq] = [save_fft]
            mfe_to_plot[len_seq] = [save_mfe]
            fftb_to_plot[len_seq] = [save_fftb]
            mx_to_plot[len_seq] = [save_mx]

list_len = list(fft_to_plot.keys())
list_len.sort()
nb_el = len(list_len)

fft_mean = [mean(fft_to_plot[len_]) for len_ in list_len]
mfe_mean = [mean(mfe_to_plot[len_]) for len_ in list_len]
fftb_mean = [mean(fftb_to_plot[len_]) for len_ in list_len]
mx_mean = [mean(mx_to_plot[len_]) for len_ in list_len]
print(len(fft_all_v))
print(len(mx_all_v))
print(mean(fft_mean))
print(ttest_ind(mx_mean, fft_mean))

Extract loop composition from known structures

from utils import paired_positions
from numpy import mean, std, median, linspace, arange, array
import numpy as np
from scipy.stats import kde, ttest_ind
from RNA import fold_compound, bp_distance
from utility.utils_analysis import get_loop_content, read_csv, read_true_struct
from math import isnan
import matplotlib.pyplot as plt
from matplotlib import rcParams, artist
from matplotlib.patches import Patch
plt.rcParams["font.family"] = "serif"
plt.rcParams["font.size"] = 15
plt.rcParams["figure.figsize"] = 9, 4

fft_pred = read_csv("./benchmark_results/fft_100n_50ms_best_nrj_scores.csv")
fftb_pred = read_csv("./benchmark_results/fft_100n_50ms_scores.csv")
vrna_pred = read_csv("./benchmark_results/mfe_scores.csv")
true_str = read_true_struct()
fft_all_v, mfe_all_v, fftb_all_v, mx_all_v = [], [], [], []
fft_to_plot, mfe_to_plot, fftb_to_plot, mx_to_plot = {}, {}, {}, {}
results = []
struct_list = []

for seq in true_str:
    if seq in fft_pred and seq in vrna_pred and seq in fftb_pred:
        len_seq, fft_st, fft_nrj, fft_nbp, fft_pvv, fft_sens = fft_pred[seq]
        len_seq, mfe_st, mfe_nrj, mfe_nbp, mfe_pvv, mfe_sens = vrna_pred[seq]
        len_seq, fftb_st, fftb_nrj, fftb_nbp, fftb_pvv, fftb_sens = fftb_pred[seq]
        save_fft, save_mfe, save_line = fft_pvv, mfe_pvv, fftb_pvv

        true_struct, name = true_str[seq]
        struct_list += [true_struct]
        int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(true_struct)
        results += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]

        fft_all_v += [save_fft]
        mfe_all_v += [save_mfe]
        fftb_all_v += [save_line]

results = array(results)
results -= results.mean(axis=0)
cov = np.cov(results.T)/results.shape[0]
V, W = np.linalg.eig(cov)
idx = V.argsort()[::-1]
W = W[:, idx]
fig = plt.figure(1)
fig.subplots_adjust(left=0.10, bottom=0.15, right=0.99, top=0.90, wspace=0.12, hspace=0.05)
# fig.subplots_adjust(left=0.06, bottom=0.13, right=0.99, top=0.90, wspace=0.12, hspace=0.1)
fig.tight_layout()

fft_f = fig.add_subplot(122)
mfe_f = fig.add_subplot(121)
# fft_f.set_aspect("equal", adjustable="box")
fft_f.grid(True, color="grey",linestyle="--", linewidth=0.2)
# mfe_f.set_aspect("equal", adjustable="box")
mfe_f.grid(True, color="grey",linestyle="--", linewidth=0.2)

wrong_fft = [(i, el) for i, el in enumerate(fft_all_v) if results.dot(W[:, 1])[i] < -0.3]
wrong_fft.sort(key=lambda el: results.dot(W[:, 1])[el[0]])
wrong_fft = [i for i, el in wrong_fft[:10]]
for el in wrong_fft:
    print(struct_list[el])
# wrong_fft = [i for i, el in enumerate(fft_all_v) if el <= 0]
wrong_mfe = [i for i, el in enumerate(mfe_all_v) if el <= 0]

proj = results.dot(W)

fft_f.scatter(results.dot(W[:,0]) , results.dot(W[:,1]) , color="grey", alpha=0.5, s=10)
fft_f.scatter(results[wrong_fft, :].dot(W[:,0]) , results[wrong_fft, :].dot(W[:,1]) , color="deepskyblue", alpha=0.5, s=14)
fft_f.arrow(0, 0, W[0,0]*0.27, W[0,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("I", xy=(W[0,0]*0.27+0.02, W[0,1]*0.27+0.02), size=15)
fft_f.arrow(0, 0, W[1,0]*0.27, W[1,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("S", xy=(W[1,0]*0.27-0.05, W[1,1]*0.27-0.02), size=15)
fft_f.arrow(0, 0, W[2,0]*0.27, W[2,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("M", xy=(W[2,0]*0.27-0.01, W[2,1]*0.27+0.02), size=15)
fft_f.arrow(0, 0, W[3,0]*0.27, W[3,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("H", xy=(W[3,0]*0.27-0.02, W[3,1]*0.27+0.03), size=15)
fft_f.arrow(0, 0, W[4,0]*0.27, W[4,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("B", xy=(W[4,0]*0.27-0.00, W[4,1]*0.27+0.03), size=15)
fft_f.arrow(0, 0, W[5,0]*0.27, W[5,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("E", xy=(W[5,0]*0.27-0.05, W[5,1]*0.27+0.00), size=15)
fft_f.set_xlabel("PC1", fontsize=17)
fft_f.tick_params(axis="y", labelleft=False, size=0)
fft_f.set_title(f"RAFFT", fontsize=20)

mfe_f.scatter(results.dot(W[:,0]) , results.dot(W[:,1]) , color="grey", alpha=0.5, s=10)
mfe_f.scatter(results[wrong_mfe, :].dot(W[:,0]) , results[wrong_mfe, :].dot(W[:,1]) , color="orangered", alpha=0.5, s=14)
mfe_f.arrow(0, 0, W[0,0]*0.27, W[0,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("I", xy=(W[0,0]*0.27+0.02, W[0,1]*0.27+0.02), size=15)
mfe_f.arrow(0, 0, W[1,0]*0.27, W[1,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("S", xy=(W[1,0]*0.27-0.05, W[1,1]*0.27-0.02), size=15)
mfe_f.arrow(0, 0, W[2,0]*0.27, W[2,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("M", xy=(W[2,0]*0.27-0.01, W[2,1]*0.27+0.02), size=15)
mfe_f.arrow(0, 0, W[3,0]*0.27, W[3,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("H", xy=(W[3,0]*0.27-0.02, W[3,1]*0.27+0.03), size=15)
mfe_f.arrow(0, 0, W[4,0]*0.27, W[4,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("B", xy=(W[4,0]*0.27-0.00, W[4,1]*0.27+0.03), size=15)
mfe_f.arrow(0, 0, W[5,0]*0.27, W[5,1]*0.27, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("E", xy=(W[5,0]*0.27-0.05, W[5,1]*0.27+0.00), size=15)
# mfe_f.tick_params(axis="y", labelleft=False, size=0)
mfe_f.set_xlabel("PC1", fontsize=17)
mfe_f.set_ylabel("PC2", fontsize=17)
mfe_f.set_title(f"MFE", fontsize=20)

# plt.savefig("img/pca_known.png", dpi=300, transparent=True)
plt.show()

from utils import paired_positions
from numpy import mean, std, median, linspace, arange, array
import numpy as np
from scipy.stats import kde, ttest_ind
from RNA import fold_compound, bp_distance
from utility.utils_analysis import get_loop_content, read_csv, read_true_struct
from math import isnan, log
import matplotlib.pyplot as plt
from matplotlib import rcParams, artist
from matplotlib.patches import Patch
plt.rcParams["font.family"] = "serif"
plt.rcParams["font.weight"] = "bold"
plt.rcParams["font.size"] = 6.7
plt.rcParams["figure.figsize"] = 7, 4.5

fft_pred = read_csv("./benchmark_results/fft_100n_50ms_best_nrj_scores.csv")
mx_pred = read_csv("./benchmark_results/mxfold_scores.csv")
fftb_pred = read_csv("./benchmark_results/fft_100n_50ms_scores.csv")
vrna_pred = read_csv("./benchmark_results/mfe_scores.csv")
true_str = read_true_struct()
fft_all_v, mfe_all_v, fftb_all_v, mx_all_v = [], [], [], []
fft_to_plot, mfe_to_plot, fftb_to_plot, mx_to_plot = {}, {}, {}, {}
results_mfe, results_fft, results_mx, results = [], [], [], []
struct_list = []

for seq in true_str:
    if seq in fft_pred and seq in vrna_pred and seq in fftb_pred and seq in mx_pred:
        len_seq, fft_st, fft_nrj, fft_nbp, fft_pvv, fft_sens = fft_pred[seq]
        len_seq, mfe_st, mfe_nrj, mfe_nbp, mfe_pvv, mfe_sens = vrna_pred[seq]
        len_seq, mx_st, mx_nrj, mx_nbp, mx_pvv, mx_sens = mx_pred[seq]
        len_seq, fftb_st, fftb_nrj, fftb_nbp, fftb_pvv, fftb_sens = fftb_pred[seq]
        save_fft, save_mfe, save_line = fft_pvv, mfe_pvv, fftb_pvv

        true_struct, name = true_str[seq]
        struct_list += [true_struct]
        int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(true_struct)
        results += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]

        int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(fft_st)
        results_fft += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]
        
        int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(mfe_st)
        results_mfe += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]
        
        int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(mx_st)
        results_mx += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]

from scipy.stats import kde
from scipy import stats
fig = plt.figure(1)
entro = lambda el_l: -sum([el * log(el) for el in el_l if el > 0.0])
# entro = lambda el_l: sum([1.0 for el in el_l if el > 0.1])
# entro = lambda el_l: el_l[0] + el_l[1]
test_cond = lambda i: len(struct_list[i]) > 80
ent_l = [entro(res) for i, res in enumerate(results) if test_cond(i)]
ent_lfft = [entro(res) for i, res in enumerate(results_fft) if test_cond(i)]
ent_lmfe = [entro(res) for i, res in enumerate(results_mfe) if test_cond(i)]
ent_lmx = [entro(res) for i, res in enumerate(results_mx) if test_cond(i)]

print(len(ent_l))
print(np.mean(ent_l))
print(np.mean(ent_lfft))
print(np.mean(ent_lmfe))

x_l = np.linspace(min(ent_l), max(ent_l), 200)

dens_ent = kde.gaussian_kde(ent_l)
dens_ent_fft = kde.gaussian_kde(ent_lfft)
dens_ent_mfe = kde.gaussian_kde(ent_lmfe)
dens_ent_mx = kde.gaussian_kde(ent_lmx)

fft_f = fig.add_subplot(111)
fft_f.plot(x_l, dens_ent(x_l), c="grey")
fft_f.plot(x_l, dens_ent_fft(x_l), c="skyblue")
fft_f.plot(x_l, dens_ent_mfe(x_l), c="orangered")
fft_f.plot(x_l, dens_ent_mx(x_l), c="blueviolet")
plt.show()

Extract loop composition from predicted structures

from utils import paired_positions
from numpy import mean, std, median, linspace, arange, array
import numpy as np
from scipy.stats import kde, ttest_ind
from RNA import fold_compound, bp_distance
from utility.utils_analysis import get_loop_content, read_csv, read_true_struct
from math import isnan
import matplotlib.pyplot as plt
from matplotlib import rcParams, artist
from matplotlib.patches import Patch
plt.rcParams["font.family"] = "serif"
plt.rcParams["font.size"] = 15
plt.rcParams["figure.figsize"] = 13, 4

fft_pred = read_csv("./benchmark_results/fft_100n_50ms_best_nrj_scores.csv")
mx_pred = read_csv("./benchmark_results/mxfold_scores.csv")
vrna_pred = read_csv("./benchmark_results/mfe_scores.csv")
true_str = read_true_struct()
fft_all_v, mfe_all_v, fftb_all_v, mx_all_v = [], [], [], []
fft_to_plot, mfe_to_plot, fftb_to_plot, mx_to_plot = {}, {}, {}, {}
results = []
struct_list = []

results_fft, results_mfe, results_mx = [], [], []
with open("./scratch/fft_loop_content.csv", "w") as out:
    out.write(f"pcc_fft,int_l,sta_l,mul_l,hai_l\n")
    for seq in true_str:
        if seq in mx_pred:
            len_seq, fft_st, fft_nrj, fft_nbp, fft_pvv, fft_sens = fft_pred[seq]
            len_seq, mfe_st, mfe_nrj, mfe_nbp, mfe_pvv, mfe_sens = vrna_pred[seq]
            len_seq, mx_st, mx_nrj, mx_nbp, mx_pvv, mx_sens = mx_pred[seq]
            true_struct, name = true_str[seq]
            int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(mfe_st)
            results_mfe += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]
            int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(fft_st)
            results_fft += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]
            int_l, sta_l, mul_l, hai_l, ext_l, bul_l = get_loop_content(mx_st)
            results_mx += [[int_l, sta_l, mul_l, hai_l, bul_l, ext_l]]
            mx_all_v += [mx_pvv]
            struct_list += [true_struct]


results_fft = array(results_fft)
results_fft -= results_fft.mean(axis=0)
cov_fft = np.cov(results_fft.T)/results_fft.shape[0]
V_fft, W_fft = np.linalg.eig(cov_fft)
idx_fft = V_fft.argsort()[::-1]
W_fft = W_fft[:, idx_fft]

results_mfe = array(results_mfe)
results_mfe -= results_mfe.mean(axis=0)
cov_mfe = np.cov(results_mfe.T)/results_mfe.shape[0]
V_mfe, W_mfe = np.linalg.eig(cov_mfe)
idx_mfe = V_mfe.argsort()[::-1]
W_mfe = W_mfe[:, idx_mfe]
# for symmetry consistency
W_mfe[:, 0] = -W_mfe[:, 0]

results_mx = array(results_mx)
results_mx -= results_mx.mean(axis=0)
cov_mx = np.cov(results_mx.T)/results_mx.shape[0]
V_mx, W_mx = np.linalg.eig(cov_mx)
idx_mx = V_mx.argsort()[::-1]
W_mx = W_mx[:, idx_mx]

fig = plt.figure(1)
fig.subplots_adjust(left=0.06, bottom=0.13, right=0.999, top=0.90, wspace=0.12, hspace=0.1)
fig.tight_layout()

fft_f = fig.add_subplot(133)
mfe_f = fig.add_subplot(132)
mx_f = fig.add_subplot(131)
fft_f.set_aspect("equal", adjustable="box")
fft_f.grid(True, color="grey",linestyle="--", linewidth=0.2)
fft_f.set_xlim([-0.4, 0.4]); fft_f.set_ylim([-0.4, 0.3])
mfe_f.set_aspect("equal", adjustable="box")
mfe_f.grid(True, color="grey",linestyle="--", linewidth=0.2)
mfe_f.set_xlim([-0.4, 0.4]); mfe_f.set_ylim([-0.4, 0.3])
mx_f.set_aspect("equal", adjustable="box")
mx_f.grid(True, color="grey",linestyle="--", linewidth=0.2)
mx_f.set_xlim([-0.4, 0.4]); mx_f.set_ylim([-0.4, 0.3])

wrong_mx = [(i, el) for i, el in enumerate(mx_all_v) if results_mx.dot(W_mx[:, 1])[i] > 0.1]
wrong_mx.sort(key=lambda el: results_mx.dot(W_mx[:, 0])[el[0]])
wrong_mx = [i for i, el in wrong_mx[:10]]
for el in wrong_mx:
    print(struct_list[el])
    
fft_f.scatter(results_fft.dot(W_fft[:,0]) , results_fft.dot(W_fft[:,1]) , color="grey", alpha=0.5, s=8)
fft_f.arrow(0, 0, W_fft[0,0]*0.3, W_fft[0,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("I", xy=(W_fft[0,0]*0.3-0.02, W_fft[0,1]*0.3+0.02), size=15)
fft_f.arrow(0, 0, W_fft[1,0]*0.3, W_fft[1,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("S", xy=(W_fft[1,0]*0.3-0.00, W_fft[1,1]*0.3+0.04), size=15)
fft_f.arrow(0, 0, W_fft[2,0]*0.3, W_fft[2,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("M", xy=(W_fft[2,0]*0.3-0.01, W_fft[2,1]*0.3+0.02), size=15)
fft_f.arrow(0, 0, W_fft[3,0]*0.3, W_fft[3,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("H", xy=(W_fft[3,0]*0.3-0.02, W_fft[3,1]*0.3+0.03), size=15)
fft_f.arrow(0, 0, W_fft[4,0]*0.3, W_fft[4,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("B", xy=(W_fft[4,0]*0.3+0.02, W_fft[4,1]*0.3+0.03), size=15)
fft_f.arrow(0, 0, W_fft[5,0]*0.3, W_fft[5,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
fft_f.annotate("E", xy=(W_fft[5,0]*0.3+0.01, W_fft[5,1]*0.3-0.03), size=15)
fft_f.tick_params(axis="y", labelleft=False, size=0)
fft_f.set_xlabel("PC1", fontsize=15)
fft_f.set_title(f"RAFFT", fontsize=20)

mfe_f.scatter(results_mfe.dot(W_mfe[:,0]) , results_mfe.dot(W_mfe[:,1]) , color="grey", alpha=0.5, s=8)
mfe_f.arrow(0, 0, W_mfe[0,0]*0.3, W_mfe[0,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("I", xy=(W_mfe[0,0]*0.3-0.02, W_mfe[0,1]*0.3+0.02), size=15)
mfe_f.arrow(0, 0, W_mfe[1,0]*0.3, W_mfe[1,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("S", xy=(W_mfe[1,0]*0.3-0.04, W_mfe[1,1]*0.3+0.03), size=15)
mfe_f.arrow(0, 0, W_mfe[2,0]*0.3, W_mfe[2,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("M", xy=(W_mfe[2,0]*0.3-0.01, W_mfe[2,1]*0.3+0.03), size=15)
mfe_f.arrow(0, 0, W_mfe[3,0]*0.3, W_mfe[3,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("H", xy=(W_mfe[3,0]*0.3-0.02, W_mfe[3,1]*0.3+0.03), size=15)
mfe_f.arrow(0, 0, W_mfe[4,0]*0.3, W_mfe[4,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("B", xy=(W_mfe[4,0]*0.3-0.06, W_mfe[4,1]*0.3+0.03), size=15)
mfe_f.arrow(0, 0, W_mfe[5,0]*0.3, W_mfe[5,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mfe_f.annotate("E", xy=(W_mfe[5,0]*0.3-0.02, W_mfe[5,1]*0.3-0.07), size=15)
mfe_f.tick_params(axis="y", labelleft=False, size=0)
mfe_f.set_xlabel("PC1", fontsize=15)
mfe_f.set_title(f"MFE", fontsize=20)

mx_f.scatter(results_mx.dot(W_mx[:,0]) , results_mx.dot(W_mx[:,1]) , color="grey", alpha=0.5, s=8)
mx_f.scatter(results_mx.dot(W_mx[:,0])[wrong_mx] , results_mx.dot(W_mx[:,1])[wrong_mx] , color="green", alpha=0.5, s=8)
mx_f.arrow(0, 0, W_mx[0,0]*0.3, W_mx[0,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mx_f.annotate("I", xy=(W_mx[0,0]*0.3-0.02, W_mx[0,1]*0.3+0.02), size=15)
mx_f.arrow(0, 0, W_mx[1,0]*0.3, W_mx[1,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mx_f.annotate("S", xy=(W_mx[1,0]*0.3+0.02, W_mx[1,1]*0.3+0.03), size=15)
mx_f.arrow(0, 0, W_mx[2,0]*0.3, W_mx[2,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mx_f.annotate("M", xy=(W_mx[2,0]*0.3-0.03, W_mx[2,1]*0.3-0.07), size=15)
mx_f.arrow(0, 0, W_mx[3,0]*0.3, W_mx[3,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mx_f.annotate("H", xy=(W_mx[3,0]*0.3-0.02, W_mx[3,1]*0.3+0.03), size=15)
mx_f.arrow(0, 0, W_mx[4,0]*0.3, W_mx[4,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mx_f.annotate("B", xy=(W_mx[4,0]*0.3+0.02, W_mx[4,1]*0.3-0.04), size=15)
mx_f.arrow(0, 0, W_mx[5,0]*0.3, W_mx[5,1]*0.3, width = 0.002, facecolor="black", head_width = 0.02)
mx_f.annotate("E", xy=(W_mx[5,0]*0.3-0.02, W_mx[5,1]*0.3+0.04), size=15)
# mx_f.tick_params(axis="y", labelleft=False, size=0)
mx_f.set_xlabel("PC1", fontsize=15)
mx_f.set_ylabel("PC2", fontsize=15)
mx_f.set_title(f"ML", fontsize=20)

# plt.savefig("img/pca_predicted.png", dpi=300, transparent=True)
plt.show()

Base pairs spanning

from utils import paired_positions
from numpy import mean, std, median, linspace, arange, array
import numpy as np
from scipy.stats import kde, ttest_ind
from RNA import fold_compound, bp_distance
from utility.utils_analysis import get_loop_content, read_csv, read_true_struct
from math import isnan
import matplotlib.pyplot as plt
from matplotlib import rcParams, artist
from matplotlib.patches import Patch
from collections import defaultdict
plt.rcParams["font.family"] = "serif"
plt.rcParams["font.weight"] = "bold"
plt.rcParams["font.size"] = 6.7
plt.rcParams["figure.figsize"] = 8, 4

fft_pred = read_csv("./benchmark_results/fft_100n_50ms_best_nrj_scores.csv")
fftb_pred = read_csv("./benchmark_results/fft_100n_50ms_scores.csv")
mx_pred = read_csv("./benchmark_results/mxfold_scores.csv")
vrna_pred = read_csv("./benchmark_results/mfe_scores.csv")
true_str = read_true_struct()
fft_all_v, mfe_all_v, fftb_all_v, mx_all_v = [], [], [], []
fft_to_plot, mfe_to_plot, fftb_to_plot, mx_to_plot = {}, {}, {}, {}
results = []
struct_list = []

results_fft = defaultdict(lambda : 0)
results_mfe = defaultdict(lambda : 0)
results_mx = defaultdict(lambda : 0)
results_norm = defaultdict(lambda : 0)
for seq, (struct, name) in true_str.items():
    seq_comp = fold_compound(seq)

    true_nb_bp = struct.count("(")
    true_nrj = seq_comp.eval_structure(struct)
    if seq in fft_pred and seq in vrna_pred and seq in mx_pred:
        len_seq, fft_st, fft_nrj, fft_nbp, fft_pvv, fft_sens = fft_pred[seq]
        len_seq, mfe_st, mfe_nrj, mfe_nbp, mfe_pvv, mfe_sens = vrna_pred[seq]
        len_seq, mx_st, mx_nrj, mx_nbp, mx_pvv, mx_sens = mx_pred[seq]
        pair_list_mfe = paired_positions(mfe_st)
        pair_list_fft = paired_positions(fft_st)
        pair_list_mx = paired_positions(mx_st)
        pair_list_true = paired_positions(struct)

        for pi, pj in pair_list_true:
            dist = abs(pi - pj)

            if (pi, pj) in pair_list_mfe:
                results_mfe[dist] += 1.0

            if (pi, pj) in pair_list_fft:
                results_fft[dist] += 1.0
                    
            if (pi, pj) in pair_list_mx:
                results_mx[dist] += 1.0

            results_norm[dist] += 1.0

dist_list = list(results_norm.keys())
dist_list.sort()
dist_list = [d for d in dist_list if results_norm[d] > 20]
# pvv_list = [mean(results[d]) for d in dist_list]
pvv_list_mfe = [results_mfe[d]/results_norm[d] for d in dist_list]
pvv_list_fft = [results_fft[d]/results_norm[d] for d in dist_list]
pvv_list_mx = [results_mx[d]/results_norm[d] for d in dist_list]

nb_el = len(dist_list)
slide_mean_mfe = [mean([d for d in pvv_list_mfe[max(0, i-10):min(nb_el, i+10)]]) for i, d_ in enumerate(pvv_list_mfe)]
slide_mean_fft = [mean([d for d in pvv_list_fft[max(0, i-10):min(nb_el, i+10)]]) for i, d_ in enumerate(pvv_list_mfe)]
slide_mean_mx = [mean([d for d in pvv_list_mx[max(0, i-10):min(nb_el, i+10)]]) for i, d_ in enumerate(pvv_list_mfe)]

fig = plt.figure(1)
fft_f = fig.add_subplot(111)

fft_f.grid(True, color="grey",linestyle="--", linewidth=0.2)
fft_f.plot(dist_list, slide_mean_mfe, linestyle="--", c="orangered", alpha = 0.5)
fft_f.scatter(dist_list, pvv_list_mfe, c="orangered", alpha = 0.8, s=6)
fft_f.plot(dist_list, slide_mean_fft, linestyle="--", c="deepskyblue", alpha = 0.5)
fft_f.scatter(dist_list, pvv_list_fft, c="deepskyblue", alpha = 0.8, s=6)
fft_f.plot(dist_list, slide_mean_mx, linestyle="--", c="blueviolet", alpha = 0.5)
fft_f.scatter(dist_list, pvv_list_mx, c="blueviolet", alpha = 0.8, s=6)
plt.xlabel("BP spanning (nb nucleotides)", fontsize=10)
plt.ylabel("% of correct BPs", fontsize=10)
# fft_f.scatter(dist_list, [results_norm[d] for d in dist_list])

fft_l = Patch(color="skyblue", label="RAFFT")
mfe_l = Patch(color="tomato", label="VRNA")
mx_l = Patch(color="blueviolet", label="ML")

name_fft = "RAFFT"
name_mfe = "MFE"
name_mx = "ML"
fft_f.legend((fft_l,mfe_l,mx_l), (f"{name_fft}", f"{name_mfe}", f"{name_mx}"), loc="lower center", ncol=4)
plt.savefig("img/bp_spanning.png")
plt.show()

Test case of folding paths

>RNA frameshift
GGGUUUGCGGUGUAAGUGCAGCCCGUCUUACACCGUGCGGCACAGGCACUAGUACUGAUGUCGUAUACAGGGCUUUUGACAU
.......((((((((..((.....)).))))))))..((.((((.(((....))))).)).))...................

seq="GGGUUUGCGGUGUAAGUGCAGCCCGUCUUACACCGUGCGGCACAGGCACUAGUACUGAUGUCGUAUACAGGGCUUUUGACAU"
python rafft.py -s $seq -n 50 -ms 5 --verbose

seq="GGGUUUGCGGUGUAAGUGCAGCCCGUCUUACACCGUGCGGCACAGGCACUAGUACUGAUGUCGUAUACAGGGCUUUUGACAU"
struct=(
    ".....(((((((((((..........)))))))))))............................................."
    "..................................((((((((................))))))))................"
    "...................(((((.............................................)))))........"
    "..................................((((.......))))................................."
    "(((((..............))))).........................................................."
)
output="img/frame_shift/step_1"
resolution=3
for i in {0..4}; do
    java -cp  ./utility/VARNAv3-93.jar fr.orsay.lri.varna.applications.VARNAcmd \
         -sequenceDBN $seq -structureDBN ${struct[$i]} -o ${output}_${i}.png -resolution $resolution\
         -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.5 -baseInner "#051C2C"\
         -baseName "#051C2C" -background "#000000" -periodNum 1000
done

seq="GGGUUUGCGGUGUAAGUGCAGCCCGUCUUACACCGUGCGGCACAGGCACUAGUACUGAUGUCGUAUACAGGGCUUUUGACAU"
struct=".......((((((((..((.....)).))))))))..((.((((.(((....))))).)).))..................."

output="img/frame_shift/wt.png"
resolution=3
java -cp  ./utility/VARNAv3-93.jar fr.orsay.lri.varna.applications.VARNAcmd \
        -sequenceDBN $seq -structureDBN $struct -o $output -resolution $resolution\
        -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.7 -baseInner "#051C2C"\
        -background "#000000" -periodNum 1000
        # -background "#000000" -periodNum 1000
# -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.5 -baseInner "#051C2C"\
# -baseName "#051C2C" -background "#000000" -periodNum 1000

Difficult structures

srp_Synt.wolf._CP000448
CCGUGCUAGAUGGGGAGGUAGCGGUGCCCUGUAACCUGCAACCCGCUAUAGCAGGGUCGAAUUCC
0.0 0.0 0.0
((.(((((....(((((((.((((....)))).))))....)))....))))).)).........
((.(((((..(((.(.(((.(((((........))).)).)))).)))))))).)).........
((.(((((..(((.(.(((.((((.(........))))).)))).)))))))).)).........
........(((((((.(.(...).).))))....(((((...........)))))))).......

srp_Heli.exil._EE659686
AAUGCGUUAGGCUGGUUUCACAGAGCUGCGAGAACCUCACGCUCUACACAGUGUAAGGAUUACA
0.0 0.0 42.86
....(.(((.((((((.....(((((.(.(((...))).)))))))).)))).))).)......
....(.(((.((((((.....(((((.(.(((...))).)))))))).)))).))).)......
....(((....(((......)))....)))....(((.(((((......))))).)))......
..(.(((....(((......)))....))).)................................

srp_Sacc.eryt._AM420293
CUGAACCCCCCCAGGGCCGGAAGGCAGCAAGGGUAGGUGGGCCCUGGCGGGUG
0.0 0.0 46.15
.......(((((((((((.(...((.......))...).)))))))).)))..
.......(((((((((((.(...((.......))...).)))))))).)))..
....((((((.....(((....))).....)))..)))..((((....)))).
((.....(((.....(((....))).....)))..))................

srp_Meth.mari._CP000745
UGGCUAGGCUGGGAAGUUAGGCGUUUCCUGUAACUCGAAAUCGCCUUUGCGAGAGCCGAAAACUUGAGGGCGGUUUUAAAUUCUGUCAUUCAUUCUCAAGUUUUGUGUAGACAUUUCGUCCUUUGGGGUAAGAUGGAGGAGGAACCUUUUUUGGAAGAAAAAGACAAACCUCCCUUAUCUUUCGAACCCCGUCAGGCCCGGAAGGGAGCAGCGGUAGA
19.4 22.86 30.3
..(((...(((((.....(((((((((.........)))).)))))..((((((((..((((((((((((.((................)).))))))))))))..)).....))))))(((.(((((((((((((.(((((...(((((((....))))))).....))))).)))))))....))))))..)))))))).....))).........
...(((.((((.......(((((((((.........)))).)))))..((((((((..(((((((((((((((.........))))).......))))))))))..)).....))))))(((.(((((((((((((.(((((...(((((((....))))))).....))))).)))))))....))))))..)))(((....)))..))))..))).
..(((.((((.....((.(((((((((.........)))).)))))..))...)))).(((((((((((((((.........))))).......))))))))))((....)).......(((.(((((((((((((.(((((...(((((((....))))))).....))))).)))))))....))))))..)))(((....)))))).........
......(((((((((((.....)))))))....((((............))))))))..............................................................................(.(((((...(((((((....))))))).....))))).)......((.....((((....(((....)))....))))..))

tmRNA_Cyan.mero._AY286123_1-236
GGGGCUGUAAGGUGUUCGACAUAAGUUGUUGUUAUUCAAUAAUUGCAAAUCAAAUUCUACCUUUUUCUAUUCCCGUUAAACAUCUAGCUGUUUAAAUUUAAAGACAAUUUAAAGACAAUCUAAAGACUCAAGAGACAAAAAUUUUAGGUAUUUAGGUAACUUUUAGGUAACUAAGUAAGGAACUAAGUAACGAAAACAACUUAUGGACGUGGGUUCAAUUCCCACCAGCUCCAAUG
23.44 23.44 33.33
(((((((............(((((((((((((((((......((((.............................(((((((......))))))).((((((.....))))))....(((((((.(((.((((........)))).))).))))))).....((((....)))))))).......))))))...)))))))))))...(((((.......))))))))))))....
(((((((............(((((((((((((((((......((((.............................(((((((......))))))).((((((.....))))))....(((((((.(((.((((........)))).))).))))))).....((((....)))))))).......))))))...)))))))))))...(((((.......))))))))))))....
(((((((............((((((((((((.....)......................................(((((((......))))))).....................................................(((((((....(((((((....))))..)))...))))))).....)))))))))))...(((((.......))))))))))))....
(((((((............((((.(.......).).......................................................................................................................................................................)))...(((((.......))))))))))))....

seq="UGGCUAGGCUGGGAAGUUAGGCGUUUCCUGUAACUCGAAAUCGCCUUUGCGAGAGCCGAAAACUUGAGGGCGGUUUUAAAUUCUGUCAUUCAUUCUCAAGUUUUGUGUAGACAUUUCGUCCUUUGGGGUAAGAUGGAGGAGGAACCUUUUUUGGAAGAAAAAGACAAACCUCCCUUAUCUUUCGAACCCCGUCAGGCCCGGAAGGGAGCAGCGGUAGA"
struct_fft="(((((((............(((((((((((((((((......((((.............................(((((((......))))))).((((((.....))))))....(((((((.(((.((((........)))).))).))))))).....((((....)))))))).......))))))...)))))))))))...(((((.......))))))))))))...."
struct_mfe="(((((((............(((((((((((((((((......((((.............................(((((((......))))))).((((((.....))))))....(((((((.(((.((((........)))).))).))))))).....((((....)))))))).......))))))...)))))))))))...(((((.......))))))))))))...."
struct_mle="(((((((............((((((((((((.....)......................................(((((((......))))))).....................................................(((((((....(((((((....))))..)))...))))))).....)))))))))))...(((((.......))))))))))))...."
struct_wts="(((((((............((((.(.......).).......................................................................................................................................................................)))...(((((.......))))))))))))...."

rna_name=tmRNA_Cyan
output_fft=img/illed_img/${rna_name}_fft.png
resolution=2
java -cp  ./utility/VARNAv3-93.jar fr.orsay.lri.varna.applications.VARNAcmd \
        -sequenceDBN $seq -structureDBN $struct_fft -o $output_fft -resolution $resolution\
        -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.7 -baseInner "#051C2C"\
        -background "#000000" -periodNum 1000

output_mfe=img/illed_img/${rna_name}_mfe.png
resolution=2
java -cp  ./utility/VARNAv3-93.jar fr.orsay.lri.varna.applications.VARNAcmd \
     -sequenceDBN $seq -structureDBN $struct_mfe -o $output_mfe -resolution $resolution\
     -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.7 -baseInner "#051C2C"\
     -background "#000000" -periodNum 1000

output_mle=img/illed_img/${rna_name}_mle.png
resolution=2
java -cp  ./utility/VARNAv3-93.jar fr.orsay.lri.varna.applications.VARNAcmd \
     -sequenceDBN $seq -structureDBN $struct_mle -o $output_mle -resolution $resolution\
     -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.7 -baseInner "#051C2C"\
     -background "#000000" -periodNum 1000

output_wts=img/illed_img/${rna_name}_wts.png
resolution=2
java -cp  ./utility/VARNAv3-93.jar fr.orsay.lri.varna.applications.VARNAcmd \
     -sequenceDBN $seq -structureDBN $struct_wts -o $output_wts -resolution $resolution\
     -algorithm naview -bpStyle "line" -fillBases True -spaceBetweenBases 0.7 -baseInner "#051C2C"\
     -background "#000000" -periodNum 1000