{"metadata":{"kernelspec":{"language":"python","display_name":"Python 3","name":"python3"},"language_info":{"name":"python","version":"3.10.12","mimetype":"text/x-python","codemirror_mode":{"name":"ipython","version":3},"pygments_lexer":"ipython3","nbconvert_exporter":"python","file_extension":".py"},"kaggle":{"accelerator":"none","dataSources":[{"sourceId":51294,"databundleVersionId":7331882,"sourceType":"competition"}],"dockerImageVersionId":30558,"isInternetEnabled":true,"language":"python","sourceType":"notebook","isGpuEnabled":false}},"nbformat_minor":4,"nbformat":4,"cells":[{"cell_type":"markdown","source":"# RNA Secondary Structure Prediction: EDA and BPPs Analysis.\n- *Arnold Bhebhe*\n- *Honors Capstone Project*\n- *Alabama State University'24*","metadata":{}},{"cell_type":"markdown","source":"### Modeule Imports","metadata":{}},{"cell_type":"code","source":"# Library imports \nimport numpy as np # linear algebra\nimport pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)\nimport matplotlib.pyplot as plt\nimport seaborn as sns\nplt.style.use('fivethirtyeight')\n\nimport warnings\nwarnings.filterwarnings(\"ignore\")","metadata":{"execution":{"iopub.status.busy":"2024-04-01T22:44:26.083076Z","iopub.execute_input":"2024-04-01T22:44:26.083686Z","iopub.status.idle":"2024-04-01T22:44:27.950637Z","shell.execute_reply.started":"2024-04-01T22:44:26.083653Z","shell.execute_reply":"2024-04-01T22:44:27.949018Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Data Loading ","metadata":{}},{"cell_type":"code","source":"# Read the CSV file 'train_data.csv' from the '/kaggle/input/stanford-ribonanza-rna-folding/' directory\n# and store it in a pandas DataFrame named 'df'\ndf = pd.read_csv('/kaggle/input/stanford-ribonanza-rna-folding/train_data.csv')\n\n# Display the first 5 rows of the DataFrame 'df'\ndf.head()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T22:44:27.952573Z","iopub.execute_input":"2024-04-01T22:44:27.953079Z","iopub.status.idle":"2024-04-01T22:46:26.281760Z","shell.execute_reply.started":"2024-04-01T22:44:27.953036Z","shell.execute_reply":"2024-04-01T22:46:26.280677Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Generate Summary Statistics for the data frame\ndf.describe().style.background_gradient(cmap='summer')","metadata":{"execution":{"iopub.status.busy":"2024-04-01T22:46:26.283383Z","iopub.execute_input":"2024-04-01T22:46:26.283964Z","iopub.status.idle":"2024-04-01T22:46:57.284957Z","shell.execute_reply.started":"2024-04-01T22:46:26.283933Z","shell.execute_reply":"2024-04-01T22:46:57.283911Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Exploratory Data Analysis (EDA)📊","metadata":{}},{"cell_type":"markdown","source":"> #### RNA Base Compistion:  distribution of different bases A, C, G, U in the sequences.","metadata":{}},{"cell_type":"code","source":"# Calculate base counts:\n# Count base occurrences in each sequence, and \n# Sum the base counts across all sequences\nbase_counts = df['sequence'].apply(lambda x: pd.Series(list(x)).value_counts()).sum()\n\n# Create a bar chart of the base counts \nbase_counts.plot(kind='bar')\n\n# Customize the plot labels and title \nplt.xlabel('Base', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Count', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.title('Base Composition', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot as an image \nplt.savefig('Base Composition.png')\n\n# Show the plot\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:29:49.397596Z","iopub.execute_input":"2024-04-01T23:29:49.398562Z","iopub.status.idle":"2024-04-01T23:45:34.394504Z","shell.execute_reply.started":"2024-04-01T23:29:49.398517Z","shell.execute_reply":"2024-04-01T23:45:34.391994Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Sequence Length Distribution: distribution of sequence lenghts of the RNAs. ","metadata":{}},{"cell_type":"code","source":"# calculate the sequence lenghts \nsequence_lengths = df['sequence'].apply(len)\n\n# Create histogram with kernel density estimation \nsns.histplot(sequence_lengths, kde=True)\nplt.xlabel('Sequence Length', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Frequency', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.title('Sequence Length Distribution', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot as an image\nplt.savefig('Sequence Length Distribution.png')\n\n# Show the plot\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:01:57.809607Z","iopub.execute_input":"2024-04-01T23:01:57.809931Z","iopub.status.idle":"2024-04-01T23:02:06.280339Z","shell.execute_reply.started":"2024-04-01T23:01:57.809903Z","shell.execute_reply":"2024-04-01T23:02:06.279205Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Experiment Type Pie Chart: distribution of different experiements","metadata":{}},{"cell_type":"code","source":"# Count occurences of each experiment type\nexperiment_type_counts = df['experiment_type'].value_counts()\n\n# Create pie chart to visualize experiement type distribution \nplt.pie(experiment_type_counts, labels=experiment_type_counts.index, autopct='%1.1f%%')\nplt.title('Experiment Type Distribution', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot\nplt.savefig('Experiment Type Distribution.png')\n\n# Show the plot\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:02:06.281583Z","iopub.execute_input":"2024-04-01T23:02:06.281911Z","iopub.status.idle":"2024-04-01T23:02:06.581685Z","shell.execute_reply.started":"2024-04-01T23:02:06.281881Z","shell.execute_reply":"2024-04-01T23:02:06.580301Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Signal-to-Noise Scatter Plot","metadata":{}},{"cell_type":"code","source":"# Determine how the number of reads affects the signal quality. \n\n# Create a scatter plot \nplt.scatter(df['signal_to_noise'], df['reads'])\n\n# Label the axes\nplt.xlabel('Signal to Noise', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Reads', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.title('Signal to Noise vs Reads', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot\nplt.savefig('Signal to Noise vs Reads.png')\n\n# Show the plot\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:02:06.583397Z","iopub.execute_input":"2024-04-01T23:02:06.583923Z","iopub.status.idle":"2024-04-01T23:02:15.857593Z","shell.execute_reply.started":"2024-04-01T23:02:06.583881Z","shell.execute_reply":"2024-04-01T23:02:15.856537Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Reactivity Heatmap: show the reactivity values across different positions","metadata":{}},{"cell_type":"code","source":"# Select only the first 50 reactivity columns\nreactivity_columns = [col for col in df.columns if col.startswith('reactivity')][:100]\nreactivity_data = df[reactivity_columns]\n\n# Convert the reactivity data to a numpy array\nreactivity_array = reactivity_data.values\n\n# Create a heatmap\nplt.figure(figsize=(10, 8))\nsns.heatmap(reactivity_array, cmap='viridis', cbar=True, xticklabels=20, yticklabels=False)\nplt.xlabel('Position', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Sample', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.title('Reactivity Heatmap (First 50 Columns)', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot\nplt.savefig('Reactivity Heatmap (First 100 Columns).png')\n\n# Show the plot\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:02:15.859114Z","iopub.execute_input":"2024-04-01T23:02:15.859437Z","iopub.status.idle":"2024-04-01T23:04:47.952131Z","shell.execute_reply.started":"2024-04-01T23:02:15.859408Z","shell.execute_reply":"2024-04-01T23:04:47.950830Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Reactivity Error Bars ","metadata":{}},{"cell_type":"code","source":"mean_reactivity = reactivity_data.mean(axis=0)\nstd_reactivity = reactivity_data.std(axis=0)\n\nx = np.arange(len(mean_reactivity))\nplt.errorbar(x, mean_reactivity, yerr=std_reactivity, fmt='o')\nplt.xlabel('Position', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Mean Reactivity', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.title('Reactivity with Error Bars', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot\nplt.savefig('Reactivity with Error Bars.png')\n\n# Show the plot\nplt.show()\n","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:47.954131Z","iopub.execute_input":"2024-04-01T23:04:47.954554Z","iopub.status.idle":"2024-04-01T23:04:50.105765Z","shell.execute_reply.started":"2024-04-01T23:04:47.954516Z","shell.execute_reply":"2024-04-01T23:04:50.104986Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"mean_reactivity = reactivity_data.mean(axis=0)\nstd_reactivity = reactivity_data.std(axis=0)\n\nx = np.arange(len(mean_reactivity))\nplt.errorbar(x, mean_reactivity, yerr=std_reactivity, fmt='o')\nplt.xlabel('Position', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Mean Reactivity', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.title('Reactivity with Error Bars', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\n\n# Save the plot\nplt.savefig('Reactivity with Error Bars.png')\n\n# Show the plot\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:50.106936Z","iopub.execute_input":"2024-04-01T23:04:50.107663Z","iopub.status.idle":"2024-04-01T23:04:52.229952Z","shell.execute_reply.started":"2024-04-01T23:04:50.107632Z","shell.execute_reply":"2024-04-01T23:04:52.229132Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Perfomance Measures for Seconadry Structure Prediction: \n- Sensitivity, Positive Predictive Value (PPV), Matthews Correlation Coefficient (MCC), and F-measure.\n- These measures collectively provide a comprehensive evaluation of a prediction model's performance.","metadata":{}},{"cell_type":"markdown","source":"> #### Sensitivity: True Positive Rate or Recall\n> - Sensitivity measures the proportion of actual positive cases that were correctly predicted as positive by a model. In RNA secondary structure prediction, it indicates how well the model identifies true base pairs.\n","metadata":{}},{"cell_type":"code","source":"def sensitivity(true_positives, false_negatives):\n    return true_positives / (true_positives + false_negatives)\n\ntrue_positives = 100  \nfalse_negatives = 20  \n\nprint(f\"True Positives: {true_positives}\")\nprint(f\"False Negatives: {false_negatives}\")\n\nsensitivity_score = sensitivity(true_positives, false_negatives)\nprint(f\"Sensitivity: {sensitivity_score}\")","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:52.231211Z","iopub.execute_input":"2024-04-01T23:04:52.231706Z","iopub.status.idle":"2024-04-01T23:04:52.237040Z","shell.execute_reply.started":"2024-04-01T23:04:52.231677Z","shell.execute_reply":"2024-04-01T23:04:52.236302Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Positive Predictive Value (PPV) or Precision:\n> - PPV measures the proportion of true positive predictions out of all positive predictions made by the model. In RNA secondary structure prediction, it indicates how accurate the positive predictions are.","metadata":{}},{"cell_type":"code","source":"def ppv(true_positives, false_positives):\n    return true_positives / (true_positives + false_positives)\n\ntrue_positives = 100  \nfalse_positives = 30  \n\nppv_score = ppv(true_positives, false_positives)\nprint(f\"Positive Predictive Value (PPV): {ppv_score}\")","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:52.238134Z","iopub.execute_input":"2024-04-01T23:04:52.238620Z","iopub.status.idle":"2024-04-01T23:04:52.254050Z","shell.execute_reply.started":"2024-04-01T23:04:52.238593Z","shell.execute_reply":"2024-04-01T23:04:52.252893Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### Matthews Correlation Coefficient (MCC)\n> - MCC takes into account true positives, true negatives, false positives, and false negatives and is particularly useful when dealing with imbalanced datasets. It ranges from -1 to 1, where 1 indicates a perfect prediction.","metadata":{}},{"cell_type":"code","source":"def mcc(true_positives, true_negatives, false_positives, false_negatives):\n    numerator = (true_positives * true_negatives) - (false_positives * false_negatives)\n    denominator = ((true_positives + false_positives) * (true_positives + false_negatives) * \n                   (true_negatives + false_positives) * (true_negatives + false_negatives)) ** 0.5\n    return numerator / denominator if denominator != 0 else 0\n\n\ntrue_positives = 100  \ntrue_negatives = 50   \nfalse_positives = 30  \nfalse_negatives = 20  \n\nmcc_score = mcc(true_positives, true_negatives, false_positives, false_negatives)\nprint(f\"Matthews Correlation Coefficient (MCC): {mcc_score}\")\n","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:52.255351Z","iopub.execute_input":"2024-04-01T23:04:52.255646Z","iopub.status.idle":"2024-04-01T23:04:52.266512Z","shell.execute_reply.started":"2024-04-01T23:04:52.255619Z","shell.execute_reply":"2024-04-01T23:04:52.265591Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"> #### F-measure (F1 Score)\n> - The F-measure is the harmonic mean of precision and recall (Sensitivity). It provides a balance between precision and recall. In RNA secondary structure prediction, it helps evaluate the overall performance of the model.","metadata":{}},{"cell_type":"code","source":"def f_measure(true_positives, false_positives, false_negatives):\n    precision = ppv(true_positives, false_positives)\n    recall = sensitivity(true_positives, false_negatives)\n    return 2 * (precision * recall) / (precision + recall) if (precision + recall) != 0 else 0\n\ntrue_positives = 100  \nfalse_positives = 30  \nfalse_negatives = 20  \n\nf_measure_score = f_measure(true_positives, false_positives, false_negatives)\nprint(f\"F-measure (F1 Score): {f_measure_score}\")","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:52.270757Z","iopub.execute_input":"2024-04-01T23:04:52.271176Z","iopub.status.idle":"2024-04-01T23:04:52.278517Z","shell.execute_reply.started":"2024-04-01T23:04:52.271147Z","shell.execute_reply":"2024-04-01T23:04:52.277225Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### ARNIE - Accessibility and Reactivity of Nucleotides for Information Extraction\n- ARNIE is a computational tool used for predicting the accessibility and reactivity of nucleotides within RNA molecules. \n- It enables researchers to analyze RNA structure and behavior by calculating base pair probabilities, providing valuable insights into RNA function and folding patterns.","metadata":{}},{"cell_type":"code","source":"%%capture\n%pip install arnie","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:04:52.279707Z","iopub.execute_input":"2024-04-01T23:04:52.280012Z","iopub.status.idle":"2024-04-01T23:05:05.478994Z","shell.execute_reply.started":"2024-04-01T23:04:52.279986Z","shell.execute_reply":"2024-04-01T23:05:05.477692Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%%capture\n%pip install draw_rna","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:05:05.481054Z","iopub.execute_input":"2024-04-01T23:05:05.481492Z","iopub.status.idle":"2024-04-01T23:05:18.369373Z","shell.execute_reply.started":"2024-04-01T23:05:05.481451Z","shell.execute_reply":"2024-04-01T23:05:18.367855Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# install eternafold\n!conda config --set auto_update_conda false\n!conda install -c bioconda eternafold --yes","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:05:18.372450Z","iopub.execute_input":"2024-04-01T23:05:18.373292Z","iopub.status.idle":"2024-04-01T23:07:30.507053Z","shell.execute_reply.started":"2024-04-01T23:05:18.373244Z","shell.execute_reply":"2024-04-01T23:07:30.505893Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"%env ETERNAFOLD_PATH=/opt/conda/bin/eternafold-bin\n%env ETERNAFOLD_PARAMETERS=/opt/conda/lib/eternafold-lib/parameters/EternaFoldParams.v1","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:30.508844Z","iopub.execute_input":"2024-04-01T23:07:30.509213Z","iopub.status.idle":"2024-04-01T23:07:30.518442Z","shell.execute_reply.started":"2024-04-01T23:07:30.509181Z","shell.execute_reply":"2024-04-01T23:07:30.517215Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"#  df is your DataFrame\nsequences = df['sequence'].tolist()\nsequences","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:30.519643Z","iopub.execute_input":"2024-04-01T23:07:30.520382Z","iopub.status.idle":"2024-04-01T23:07:30.796142Z","shell.execute_reply.started":"2024-04-01T23:07:30.520336Z","shell.execute_reply":"2024-04-01T23:07:30.795356Z"},"collapsed":true,"jupyter":{"outputs_hidden":true},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from arnie.mfe import mfe\nsequence =\"GGGAACGACUCGAGUAGAGUCGAAAACGUACGUGGAGACACGUACGACGAGACUUCGGUCUCGAAUAGCUCGAACCGGUGCCGAGCGCGCACGAGGCGUCGGUGGUCGGCGAGCCCACGGACGCAAAAACUCGUGCGUAAAUAAAUAGCAUUGGAAGUGACGUCGACCGUUCGCGGUCGACGUCACAAAAGAAACAACAACAACAAC\"\nstructure = mfe(sequence,package=\"eternafold\")\nprint(structure)","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:30.797706Z","iopub.execute_input":"2024-04-01T23:07:30.798005Z","iopub.status.idle":"2024-04-01T23:07:30.955622Z","shell.execute_reply.started":"2024-04-01T23:07:30.797979Z","shell.execute_reply":"2024-04-01T23:07:30.954473Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from draw_rna.ipynb_draw import draw_struct\ndraw_struct(sequence, structure)","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:30.956650Z","iopub.execute_input":"2024-04-01T23:07:30.957185Z","iopub.status.idle":"2024-04-01T23:07:33.645710Z","shell.execute_reply.started":"2024-04-01T23:07:30.957147Z","shell.execute_reply":"2024-04-01T23:07:33.644881Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from arnie.mfe import mfe\nsequence = \"GGGAACGACTCGAGTAGAGTCGAAAAATAAGAGTGATTGGCGTCCGTACGTACCCTTTCTACTCTCAAACTCTTGTTAGTTTAAATCTAATCTAAACTTTATAAACGGCACTTCCTGTGTGTCCATGCCCGTGGGCTTGGTCTTGTCATAGTGCTGACATTTGTGGTTCCTTGGTTTTTGTTCTCTGCCAGTGACGTGTCCATTCGGCGCCAGCAGCCCACCCATAGGTTGCATAATGGCAAAGATGGGCAAATACGGTCTCGGCTTCAAATGGGCCCCAGAATTTCCATGGATGCTTCCGAACGCATCGGAGAAGTTGGGTAGCCCTGAGAGGTCAGAGGAGGATGGGTTTTGCCCCTCTGCTGCGCAAGAACCAAAAACTAAAGGAAAAACTTTGATTAATCACGTGAGGGTGGGGATCCTGTTCGCAGGATCCAAAAGAAACAACAACAACAAC\"\nstructure = mfe(sequence,package=\"eternafold\")\nprint(structure)","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:33.647183Z","iopub.execute_input":"2024-04-01T23:07:33.647744Z","iopub.status.idle":"2024-04-01T23:07:34.829914Z","shell.execute_reply.started":"2024-04-01T23:07:33.647711Z","shell.execute_reply":"2024-04-01T23:07:34.828799Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from draw_rna.ipynb_draw import draw_struct\ndraw_struct(sequence, structure)","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:34.831500Z","iopub.execute_input":"2024-04-01T23:07:34.832252Z","iopub.status.idle":"2024-04-01T23:07:40.581014Z","shell.execute_reply.started":"2024-04-01T23:07:34.832209Z","shell.execute_reply":"2024-04-01T23:07:40.580018Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"from arnie.bpps import bpps\nbpps(sequence,package=\"eternafold\")","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:40.582421Z","iopub.execute_input":"2024-04-01T23:07:40.582754Z","iopub.status.idle":"2024-04-01T23:07:41.782358Z","shell.execute_reply.started":"2024-04-01T23:07:40.582725Z","shell.execute_reply":"2024-04-01T23:07:41.781223Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Example of saving BPPs to a CSV file\nimport pandas as pd\n\nbpps_data = bpps(sequence, package=\"eternafold\")\nbpps_df = pd.DataFrame(bpps_data)\nbpps_df.to_csv('bpps_data.csv')","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:41.783869Z","iopub.execute_input":"2024-04-01T23:07:41.784592Z","iopub.status.idle":"2024-04-01T23:07:43.102558Z","shell.execute_reply.started":"2024-04-01T23:07:41.784551Z","shell.execute_reply":"2024-04-01T23:07:43.101291Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"\"\"\"\nThis script calculates base pair probabilities (BPPs) for a given RNA sequence.\n\"\"\"\n\nfrom arnie.bpps import bpps\n\ndef calculate_bpps(sequence):\n    \"\"\"\n    Calculate BPPs for a given RNA sequence.\n    \n    Args:\n        sequence (str): The RNA sequence.\n        \n    Returns:\n        dict: A dictionary containing BPPs.\n    \"\"\"\n    return bpps(sequence, package=\"eternafold\")\n\n# Usage\nbpps_data = calculate_bpps(sequence)\n\nbpps_data\n","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:43.104146Z","iopub.execute_input":"2024-04-01T23:07:43.104571Z","iopub.status.idle":"2024-04-01T23:07:44.299965Z","shell.execute_reply.started":"2024-04-01T23:07:43.104530Z","shell.execute_reply":"2024-04-01T23:07:44.298884Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"# Drop columns with NaN values (if any)\nbpps_df_cleaned = bpps_df.dropna(axis=1, how='any')\n\n# Calculate Mean BPPs for Each Position\nmean_bpps = bpps_df_cleaned.mean()\n\n# Plot Mean BPPs\nplt.figure(figsize=(10, 6))\nplt.plot(mean_bpps)\nplt.title('Mean Base Pair Probabilities', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.xlabel('Position', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Mean BPP', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.savefig('Mean Base Pair Probabilities.png')\nplt.show()\n\n# Calculate Total BPP for Each Sequence\ntotal_bpps = bpps_df_cleaned.sum(axis=1)\n\n# 4. Plot Total BPP Distribution\nplt.figure(figsize=(10, 6))\nplt.hist(total_bpps, bins=30, color='skyblue', edgecolor='black')\nplt.title('Total Base Pair Probabilities Distribution', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.xlabel('Total BPP', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Frequency', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.savefig('Total Base Pair Probabilities Distribution.png')\nplt.show()\n\n\n# Visualize BPPs Heatmap for a Few Sequences\nsample_sequences = bpps_df.sample(n=5, random_state=42)\n\nplt.figure(figsize=(10, 6))\nplt.imshow(sample_sequences.values, cmap='viridis', aspect='auto')\nplt.colorbar(label='BPP')\nplt.title('Base Pair Probabilities Heatmap', fontsize = 14, fontweight = 'bold', color = 'darkgreen')\nplt.xlabel('Position', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.ylabel('Sample Sequence', fontsize = 12, fontweight = 'bold', color = 'darkblue')\nplt.savefig('Base Pair Probabilities Heatmap.png')\nplt.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:44.301444Z","iopub.execute_input":"2024-04-01T23:07:44.301774Z","iopub.status.idle":"2024-04-01T23:07:45.786077Z","shell.execute_reply.started":"2024-04-01T23:07:44.301743Z","shell.execute_reply":"2024-04-01T23:07:45.784984Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"code","source":"import plotly.graph_objects as go\n\n# Define the performance measures\nlabels = ['Sensitivity', 'Positive Predictive Value (PPV)', 'Matthews Correlation Coefficient (MCC)', 'F-measure (F1 Score)']\nscores = [sensitivity_score, ppv_score, mcc_score, f_measure_score]\n\n# Create a horizontal bar plot\nfig = go.Figure(data=[go.Bar(\n    y=labels,\n    x=scores,\n    orientation='h',\n    marker=dict(color=['blue', 'green', 'red', 'purple'])\n)])\n\nfig.update_layout(\n    title='Performance Measures',\n    xaxis_title='Score',\n    yaxis_title='Metric',\n    yaxis=dict(autorange='reversed')\n)\n\nfig.show()","metadata":{"execution":{"iopub.status.busy":"2024-04-01T23:07:45.787442Z","iopub.execute_input":"2024-04-01T23:07:45.787764Z","iopub.status.idle":"2024-04-01T23:07:46.286217Z","shell.execute_reply.started":"2024-04-01T23:07:45.787735Z","shell.execute_reply":"2024-04-01T23:07:46.285201Z"},"trusted":true},"execution_count":null,"outputs":[]},{"cell_type":"markdown","source":"### Conclusions","metadata":{}},{"cell_type":"markdown","source":"> The model demonstrates promising performance in predicting positive cases, with a sensitivity (true positive rate) of 83.3%, indicating its ability to correctly identify true positive cases. The positive predictive value (PPV) or precision score of 76.9% means that out of all the positive predictions made by the model, 76.9% were accurate.\n\n> The Matthews Correlation Coefficient (MCC) score of 47.1% suggests a moderate level of agreement between the predicted and actual classes, indicating room for improvement to achieve a higher level of agreement.\n\n> The F1 score of 80.0% reflects a good balance between precision and recall, demonstrating that the model performs well in achieving a balance between accurate positive predictions and capturing actual positives.\n\n> Overall, the model shows strong results in terms of sensitivity and F1 score, but there is potential for further refinement, especially in enhancing the overall agreement between predicted and actual classes, as indicated by the MCC score. Fine-tuning the model or exploring additional features may lead to improved predictions.","metadata":{}},{"cell_type":"markdown","source":"#### ","metadata":{}}]}