Overview

---

Modeling played a central role in the siREN development process, guiding design choices and reducing experimental uncertainty, cost and time. Our challenge required selecting appropriate siRNAs, ensuring their thermodynamic and structural compatibility with the RNAi machinery, and developing a stable, targeted lipid nanoparticle (LNP) delivery system. To address this, we developed and employed complementary computational tools across multiple scales.

1. Machine learning for efficacy prediction: We developed and trained a machine learning model on our harmonized dataset (siRBench) to predict candidate siRNA silencing efficacy before experimental validation.
2. Thermodynamic stability analysis: We performed thermodynamic analyses using RNAfold Web Server, DuplexFold, and HDOCK to predict siRNA folding and hybridization ΔG values and estimate the stability of siRNA-Ago2 binding.
3. Structural modeling: We used AlphaFold to predict the Ago2–siRNA–mRNA ternary complex structure, providing mechanistic insights into strand positioning and binding interactions.
4. Delivery system optimization: We modeled lipid nanoparticle components using Avogadro2, analyzed membrane stability with CHARMM-GUI, and employed ROSETTA InterfaceAnalyzer to select antibodies for ROR1 targeting.

The combination of these models enabled us to understand and optimize crucial aspects of siREN: from predicting the best-performing siRNAs, to designing their delivery within a stable, functionalized nanoparticle. Through a combination of data-driven predictions and molecular simulations, our modeling not only supports experimental design but also provides a framework for RNAi-based therapeutic development.

Click on each button to review our Modeling Results!

siRNA Modeling

---

Our team decided to experimentally evaluate the silencing efficacy of both previously validated siRNAs and newly designed candidates targeting the BCL2 and BTK genes. The Wet Lab team conducted literature review to identify high-performing siRNAs from published studies, while the Dry Lab team designed novel sequences. After compiling these candidates, we faced a critical challenge: selecting the most promising siRNAs for experimental validation, a task requiring integrated modeling approaches.

1

Candidate Collection

Wet Lab: literature review to identify validated siRNAs for BCL2 and BTK
Dry Lab: design of new siRNAs. Combined list prepared for modeling-based selection.

2

Machine Learning Ranking

To identify the optimal siRNAs against BCL2 and BTK, we combined machine learning predictions with thermodynamic and structural modeling. We first employed our machine learning models, trained on the harmonized siRBench dataset, to rank candidate siRNAs by predicted silencing efficacy.

3

Thermodynamic & Structural Screening

We refined these predictions using RNAfold Web Server, DuplexFold, and HDOCK to evaluate siRNA folding energies, hybridization thermodynamics, and siRNA-Ago2 binding affinity.

4

Complex Modeling (AlphaFold)

Finally, we used AlphaFold to model the Ago2–siRNA–mRNA ternary complex for top candidates, verifying proper guide strand positioning within the Ago2 binding pocket and confirming target mRNA accessibility for cleavage.

5

Handoff to Wet Lab

This multi-scale approach enabled us to predict not only which siRNAs had high silencing potential, but also their functional viability within the cellular environment. By integrating efficacy predictions with thermodynamic stability and structural compatibility assessments, we reduced experimental trial-and-error and guided our proof-of-concept studies. The modeling results were provided to the Wet Lab team for final siRNA sequence selection and experimental validation.

To make our documentation clearer and avoid repeatedly displaying full nucleotide sequences, we assigned a unique reference code to each siRNA candidate. The codes follow a simple format:

E

Candidate siRNA E1 (#Target_Gene) → for siRNAs that have been previously experimentally validated.

D

Candidate siRNA D1 (#Target_Gene) → for siRNAs that we designed ourselves for siREN.

This methodology allowed us to easily refer to each siRNA throughout our modeling and experimental workflow, while clearly distinguishing between validated and newly designed candidates.

List of All BCL2-Targeting siRNAs

Assigned Reference Code	siRNA and Exact Target mRNA Sequences
Candidate siRNA E1 (BCL2)	Antisense (5′→3′): GGUCUGCAGCGGCGAGGUCCUGGC Target (5′→3′): GCCAGGACCUCGCCGCUGCGACC
Candidate siRNA E2 (BCL2)	Antisense (5′→3′): GGCCGUACAGUUCCACAAAGGCAU Target (5′→3′): AUGCCUUUGUGGAACUGUACGGCC
Candidate siRNA E3 (BCL2)	Antisense (5′→3′): GUCUGCAGCGGCGAGGUCCUGGC Target (5′→3′): GCCAGGACCUCGCCGCUGCAGAC
Candidate siRNA E4 (BCL2)	Antisense (5′→3′): GCCGUACAGUUCCACAAAGGCAU Target (5′→3′): AUGCCUUUGUGGAACUGUACGGC
Candidate siRNA E5 (BCL2)	Antisense (5′→3′): UCAGGUACUCAGUCAUCCACAGG Target (5′→3′): CCUGUGGAUGACUGAGUACCUGA
Candidate siRNA E6 (BCL2)	Antisense (5′→3′): CGCAGCTTATAATGGATGTACT Target (5′→3′): AGUACAUCCAUUAUAAGCUGCG
Candidate siRNA E7 (BCL2)	Antisense (5′→3′): UACAGUUCCACAAAGGCAUCCUU Target (5′→3′): AAGGAUGCCUUUGUGGAACUGUA
Candidate siRNA E8 (BCL2)	Antisense (5′→3′): GAAGGGCGUCAGGUGCAGC Target (5′→3′): GAAGGGCGUCAGGUGCAGC
Candidate siRNA E9 (BCL2)	Antisense (5′→3′): GCAGGCAUGUUGACUUCAC Target (5′→3′): GUGAAGUCAACAUGCCUGC
Candidate siRNA E10 (BCL2)	Antisense (5′→3′): AAAGGCAUCCCAGCCUCCG Target (5′→3′): CGGAGGCUGGGAUGCCUUU
Candidate siRNA E11 (BCL2)	Antisense (5′→3′): AGCUUAUAAUGGAUGUACU Target (5′→3′): AGUACAUCCAUUAUAAGCU
Candidate siRNA D1 (BCL2)	Antisense (5′→3′): UCAAAGAAGGCCACAAUCCUC Target (5′→3′): GAGGAUUGUGGCCUUCUUUGAGU

List of All BTK-Targeting siRNAs

Assigned Reference Code	siRNA and Exact Target mRNA Sequences
Candidate siRNA E1 (BTK)	Antisense (5′→3′): CUCCUUGAUCGUUUACCAAUU Target (5′→3′): TTGGTAAACGATCAAGGAG
Candidate siRNA E2 (BTK)	Antisense (5′→3′): UGAAACCUCCUUCUUUCCCUU Target (5′→3′): GGGAAAGAAGGAGGTTTCA
Candidate siRNA E3 (BTK)	Antisense (5′→3′): UUCAAAGUCAUACUCAUAGUU Target (5′→3′): ACTATGAGTATGACTTTGAA
Candidate siRNA D1 (BTK)	Antisense (5′→3′): UUCAAAGUCAUACUCAUAGUA Target (5′→3′): TACTATGAGTATGACTTTGAACG
Candidate siRNA D2 (BTK)	Antisense (5′→3′): UAUUGAACCCUUCUUACUGCC Target (5′→3′): GGCAGTAAGAAGGGTTCAATAGA
Candidate siRNA D3 (BTK)	Antisense (5′→3′): UCUAUUGAACCCUUCUUACUG Target (5′→3′): CAGTAAGAAGGGTTCAATAGATG
Candidate siRNA D4 (BTK)	Antisense (5′→3′): UUUCAAUGAUUGAAAUUUGCU Target (5′→3′): AGCAAATTTCAATCATTGAAAGG
Candidate siRNA D5 (BTK)	Antisense (5′→3′): AUAGUUACUAGGAAUGUAGCC Target (5′→3′): GGCTACATTCCTAGTAACTATGT
Candidate siRNA D6 (BTK)	Antisense (5′→3′): UUUGGAAUACCACUCAUACAU Target (5′→3′): ATGTATGAGTGGTATTCCAAACA
Candidate siRNA D7 (BTK)	Antisense (5′→3′): AAUGAAACCUCCUUCUUUCCC Target (5′→3′): GGGAAAGAAGGAGGTTTCATTGT
Candidate siRNA D8 (BTK)	Antisense (5′→3′): AUAUUUGAGCCUGGAUAUGAG Target (5′→3′): CTCATATCCAGGCTCAAATATCC
Candidate siRNA D9 (BTK)	Antisense (5′→3′): AGAUUCAUCAUGACUUUGGCU Target (5′→3′): AGCCAAAGTCATGATGAATCTTT
Candidate siRNA D10 (BTK)	Antisense (5′→3′): UUACCAAACAGUUUCGAGCUG Target (5′→3′): CAGCTCGAAACTGTTTGGTAAAC
Candidate siRNA D11 (BTK)	Antisense (5′→3′): UAGUAAAUCUCUCAUAUGGCA Target (5′→3′): TGCCATATGAGAGATTTACTAAC

Machine learning is transforming synthetic biology by enabling predictive design to replace traditional trial-and-error approaches. Advanced machine learning models can generate accurate predictions that accelerate laboratory workflows, advance research, and reduce experimental costs. In siREN, our machine learning model provided initial predictions of candidate siRNA silencing efficacy and shaped our project design. Our code for feature extraction and model development is available on our GitLab Repository.

To develop our siRNA efficacy prediction model, we made the following assumptions:

We developed machine learning models to predict siRNA silencing efficiency by analyzing sequence and thermodynamic features of the siRNA antisense strand and the extended target mRNA site. The model outputs a value from 0 to 1, where 0 indicates no silencing effect and 1 indicates maximum silencing effect. Higher outputs indicate more promising candidate siRNAs. This output represents the percentage of target mRNA knockdown (0-100% corresponding to 0-1 output values).

As Dr. Andigoni Malousi emphasized during our consultation (see Integrated Human Practices Page), a rich, high-quality training dataset is essential for developing a high-performance machine learning model that produces reliable results. We searched for publicly available datasets containing siRNA antisense strands, target mRNA sites (either extended or exact site) and silencing efficiency values. Following Dr. Malousi's guidance on combining datasets, we carefully verified that all efficiency values represented silencing efficacy on the same scale. We therefore preprocessed all siRNAs to ensure they met specific standardization criteria before integration into our final dataset (siRBench).

Why a Harmonized Dataset Was Needed

Creating a high-quality dataset was an essential step before moving on to develop our machine learning models. The are multiple public siRNA sequences datasets, however when these datasets are compared to each other, some issues arise:

• Inconsistent efficacy reporting: Silencing efficiencies are experimentally measured using different assays, protocols, and cell lines, and are often reported on different scales.
• Heterogeneous sequence formats: siRNAs vary in length (commonly 19 nt or 21 nt) and in orientation (5’-3’ or 3’-5’), making direct integration challenging.
• Data inconsistencies: Several datasets contain minor errors that undermine reproducibility.

Our Goal: siRBench

As iGEM Thessaloniki 2025, we set the goal of creating siRBench, the first harmonized training dataset that comprises as many siRNAs as possible. The siRBench dataset is available on our GitLab Repository.

Dataset Composition

The final siRBench dataset comprises 4098 siRNA antisense sequences from 11 different source publications: Huesken¹, Takayuki², Shabalina³, Simone⁴,Amarzguioui⁵, Harborth ⁶, Hsieh ⁷, Khvorova ⁸, Reynolds ⁹, Vickers ¹⁰, Ui-Tei ¹¹. Note that Mixset is an aggregated dataset containing sequences from seven of these sources (Amarzguioui, Harborth, Hsieh, Khvorova, Reynolds, Vickers, and Ui-Tei), which we separated back into their original publications for proper cell line annotation.

The dropdown buttons below include the documentation of the preprocessing of each different dataset we found.

The above preprocessing procedure ensured that we created a balanced, normalized and high-quality training dataset for our model, that contains as many siRNA sequences as possible. We have uploaded our siRBench dataset on our GitLab Repository and future iGEM teams can deploy it for their own modeling endeavors!

Feature engineering transforms raw data into numerical representations that machine learning models can process. For siRBench, we converted each siRNA-target mRNA pair into a feature set capturing sequence composition, structural properties, and thermodynamic characteristics.

Machine learning algorithms cannot directly interpret nucleotide sequences, which are categorical rather than numerical data. Instead, models learn from quantitative features such as nucleotide frequencies, position-specific k-mers, predicted secondary structures, duplex stability (ΔG), and other sequence-derived properties. By pre-calculating these features, we enabled our models to indentify motifs that distinguish highly effective siRNAs from less effective ones.

Feature engineering was particularly critical for integrating heterogeneous data from multiple public datasets. It provided a unified numerical framework for comparing siRNAs measured under different experimental conditions and cell lines. The quality and diversity of the engineered features directly impact model accuracy and generalizability, making feature engineering a central pillar of reliable siRNA efficacy prediction.

Through extensive literature review, we identified 100 thermodynamic and structural features to enrich our training dataset. All feature values are rounded to three decimal places. The dropdown button below lists these features with descriptions of their biological significance.

Thermodynamics

• ends:Terminal asymmetry metric comparing 5′ and 3′ end stability of the siRNA duplex.
• DG_total: Total Gibbs free energy (ΔG) of duplex formation calculated from nearest-neighbor parameters ¹².
• DH_total: Total enthalpy change (ΔH, kcal/mol) of duplex formation from nearest-neighbor parameters ¹².
• DG_pos1..18: Per-step stacking ΔG along the duplex at each of the 18 nearest-neighbor positions in the 19-bp region.
• DH_pos1..18: Per-step stacking ΔH at corresponding positions.

Constraint energies (RNAfold/RNAcofold)

• single_energy_total: Minimum free energy (MFE) of the isolated guide strand (RNAfold).
• single_energy_pos1..19: Position-specific accessibility of each nucleotide in the guide strand (RNAfold).
• duplex_energy_total: Hybridization ΔG for the guide–target pair (RNAcofold).
• duplex_energy_sirna_pos1..19: Per-nucleotide contribution of the guide to the duplex energy at each aligned position (RNAcofold).
• duplex_energy_target_pos1..19: Per-nucleotide contribution of the target site to the duplex energy (RNAcofold).

RNAup

• RNAup_open_dG: Energy cost to unpair interacting regions on both guide and target strands prior to binding.
• RNAup_interaction_dG: Hybridization energy gained upon binding of unpaired regions.

With siRBench finalized, we developed machine learning models to predict siRNA silencing efficacy. Our model development followed iterative engineering cycles using the Design-Build-Test-Learn framework (see our Engineering Page).

We explored multiple machine learning approaches and ultimately implemented two models with distinct architectures, both fully documented in our GitLab Repository.

We sought an advanced machine learning architecture suited to our task's characteristics: multiple engineered features, a medium-sized dataset, and potentially complex non-linear relationships. After literature review, we selected LightGBM (Light Gradient Boosting Machine) ¹³, an efficient implementation of gradient boosting decision trees.

Unlike linear models, LightGBM captures non-linear feature interactions, where silencing depends on the complex interplay of thermodynamic stability, accessibility, and positional effects rather than simple linear rules. LightGBM builds an ensemble of decision trees sequentially, with each tree correcting errors from previous ones. Its leaf-wise growth strategy expands branches that maximize error reduction, enabling the model to learn highly specific patterns. Efficiency optimizations including histogram-based feature binning, gradient-based sampling, and exclusive feature bundling make it faster and more memory-efficient than other gradient boosting implementations.

Advantages for our dataset:

• Efficiently handles medium-sized datasets (~4,100 siRNAs) without requiring the massive training sets needed for deep learning.
• Naturally processes correlated and position-specific features (DG_pos1..18, DH_pos1..18, duplex energies) to identify non-linear patterns that simpler models cannot capture.

Full documentation for the LightGBM model can be found in our GitLab Repository.

After implementing and evaluating LightGBM, we expanded our modeling strategy by exploring a fundamentally different machine learning approach. Our goal was to determine whether an alternative learning mechanism could capture additional patterns and potentially improve predictive performance.

Given our tabular dataset with complex non-linear dependencies, we focused our literature review on models specifically designed for tabular data prediction. We selected the TabPFN (Tabular Prior-Data Fitted Network) regressor ¹⁴, a transformer-based model pre-trained on millions of synthetic tasks, that enables accurate predictions with minimal hyperparameter tuning.

TabPFN leverages transformer architecture and extensive pre-training to capture complex, non-linear feature interactions in tabular data without requiring dataset-specific training. This fundamentally different mechanism from gradient boosting allowed us to explore whether alternative model architectures could enhance siRNA efficacy prediction.

Full documentation for TabPFN model development is available in our GitLab Repository.

We evaluated performance using Mean Squared Error (MSE) and Pearson correlation coefficient (R) to measure alignment between predicted and actual silencing efficacies.

	LightGBM Model Score	TabPFN Regressor Score
Mean Squared Error (MSE)	0,0529	0,0665
Pearson Coefficient (R)	0,5454	0,5576

From the evaluation of the TabPFN regressor, we obtained a Mean Squared Error (MSE) of 0,0665 and a Pearson correlation coefficient (R) of 0,5576 on the test set. Although this MSE value is slightly higher than the one achieved by LightGBM (0,0529), it still reflects a reasonable prediction error. On the other hand, the Pearson correlation coefficient is slightly higher compared to LightGBM (0,5454), showing that TabPFN was able to capture the overall trend between predicted and true efficacies slightly better.

• LightGBM achieved lower prediction error, meaning its predictions were closer to the experimental values on average.
• TabPFN showed a marginally stronger correlation with the true values, suggesting it better preserved the global ranking and direction of siRNA efficacy trends.

Considering these complementary strengths, we propose both models as part of our modeling strategy, since each offers distinct advantages for understanding and predicting siRNA silencing efficacy.

The tables below show the silencing efficacy score by both proposed modes for all siRNA sequences. All silencing scores were rounded up to the 3rd decimal point. Below, we highlight the main observations and conclusions we made by studying the silencing scores of both tables:

Assigned Reference Code	LightGBM Model Score	TabPFN Regressor Score
Candidate siRNA E1 (BCL2)	0,382	0,492
Candidate siRNA E2 (BCL2)	0,391	0,481
Candidate siRNA E3 (BCL2)	0,538	0,541
Candidate siRNA E4 (BCL2)	0,509	0,528
Candidate siRNA E5 (BCL2)	0,653	0,594
Candidate siRNA E6 (BCL2)	0,633	0,583
Candidate siRNA E7 (BCL2)	0,722	0,607
Candidate siRNA E8 (BCL2)	0,487	0,530
Candidate siRNA E9 (BCL2)	0,496	0,548
Candidate siRNA E10 (BCL2)	0,543	0,529
Candidate siRNA E11 (BCL2)	0,609	0,579
Candidate siRNA D1 (BCL2)	0,753	0,611

Assigned Reference Code	LightGBM Model Score	TabPFN Regressor Score
Candidate siRNA E1 (BTK)	0,413	0,525
Candidate siRNA E2 (BTK)	0,701	0,609
Candidate siRNA E3 (BTK)	0,691	0,613
Candidate siRNA D1 (BTK)	0,691	0,613
Candidate siRNA D2 (BTK)	0,699	0,612
Candidate siRNA D3 (BTK)	0,673	0,605
Candidate siRNA D4 (BTK)	0,659	0,602
Candidate siRNA D5 (BTK)	0,699	0,612
Candidate siRNA D6 (BTK)	0,735	0,613
Candidate siRNA D7 (BTK)	0,686	0,619
Candidate siRNA D8 (BTK)	0,611	0,602
Candidate siRNA D9 (BTK)	0,688	0,609
Candidate siRNA D10 (BTK)	0,792	0,621
Candidate siRNA D11 (BTK)	0,696	0,618

BCL2-targeting siRNAs (Table 1)

• The LightGBM model predictions ranged from 0,382 (E1) to 0,753 (D1).
• The TabPFN model predictions ranged from 0,481 (E2) to 0,611 (D1).
• Both models consistently identified Candidate siRNA D1 (BCL2) as the top-performing sequence, with the highest predicted silencing efficacy.
• Among previously validated siRNAs, E7 (0,722 / 0,607) and E5 (0,653 / 0,594) ranked as the strongest candidates.
• There was good agreement between the models, with only small differences in ranking order.

BTK-targeting siRNAs (Table 2)

• The LightGBM model predictions ranged from 0,413 (E1) to 0,792 (D10).
• The TabPFN predictions ranged from 0,525 (E1) to 0,621 (D10).
• Both models identified Candidate siRNA D10 (BTK) as the most effective sequence, achieving the highest scores in both frameworks.
• Several designed siRNAs (D6, D2, D5) also scored highly, generally outperforming the previously validated E-series candidates.
• Agreement between LightGBM and TabPFN was again strong, with the same top candidates consistently emerging.

siRNA potency depends not only on sequence but also on thermodynamic properties that govern its pathway and activity within RISC. For successful gene silencing, the antisense strand must form stable, specific interactions with both the sense strand and Argonaute 2 (Ago2). However, excessive stability can inhibit strand unwinding and reduce silencing efficiency.

Thermodynamic modeling enables prediction of these interactions before experimentation. By calculating free energies using computational tools, we can assess whether a candidate siRNA will successfully enter RISC, navigate the pathway, and engage its target mRNA. This computational screening is critical for prioritizing sequences with the highest probability of successful gene silencing.

In siREN, thermodynamic modeling played a key role in selecting siRNAs with optimal stability and accessibility, ensuring only the most promising candidates advanced to structural and experimental validation.

We made the following assumption while implementing the thermodynamic modeling for siREN:

Through extensive literature review, we identified three critical thermodynamic parameters for siRNA selection:

• GC Content: Percentage of nucleotides being G or C
• Folding ΔG:Quantifies the propensity of the guide strand to form intramolecular secondary structures. Less structured strands (ΔG ≈ 0 or positive) tend to load into RISC and function more effectively, as Ago2 does not need to unwind internal hairpins before target search ¹⁵.
• Hybridization ΔG: Quantifies the free energy when the guide strand forms a duplex with its target mRNA. More negative ΔG values indicate greater duplex stability and typically correlate with improved silencing efficacy. ¹⁵
• Binding antisense strand–Ago2 ΔG: Reflects the interaction energy between the guide strand and Argonaute-2. More negative values indicate stronger, more stable binding that supports effective RISC assembly. However, binding must be strong enough to anchor the guide strand while maintaining sufficient flexibility at the PAZ domain to enable target cleavage. ¹⁵

The tables below present predicted energies for all candidate siRNAs, ranked by their thermodynamic favorability.

Assigned Reference Code	GC Content (%)	Folding DG (kcal/mol)	Hybridization DG (kcal/mol)	HDOCK Score (kcal/mol)
Candidate siRNA E5 (BCL2)	52	0	-43,2	-409,96
Candidate siRNA E4 (BCL2)	52	-3,8	-42,9	-397,27
Candidate siRNA D1 (BCL2)	48	-0,1	-36,4	-373,94
Candidate siRNA E6 (BCL2)	41	-1,3	-35,7	-365,26
Candidate siRNA E7 (BCL2)	43	-0,7	-40,1	-363,64
Candidate siRNA E9 (BCL2)	53	0	-32,9	-331,91
Candidate siRNA E11 (BCL2)	32	0	-27,3	-351,97
Candidate siRNA E1 (BCL2)	71	-4,3	-49,2	-320,77
Candidate siRNA E2 (BCL2)	54	-3,5	-46,2	-500,91
Candidate siRNA E10 (BCL2)	63	-5,1	-38,4	-355,78
Candidate siRNA E3 (BCL2)	74	-3,7	-51,5	-272,8
Candidate siRNA E8 (BCL2)	68	-0,6	-9,1	-352,79

Assigned Reference Code	GC Content (%)	Folding DG (kcal/mol)	Hybridization DG (kcal/mol)	HDOCK Score (kcal/mol)
Candidate siRNA E1 (BTK)	38	0	-31	-365,6
Candidate siRNA O2 (BTK)	43	0	-34,2	-353,57
Candidate siRNA O3 (BTK)	38	0	-32,3	-340,35
Candidate siRNA O7 (BTK)	43	0	-33,7	-329,12
Candidate siRNA O9 (BTK)	38	0	-33,6	-354,89
Candidate siRNA O11 (BTK)	33	0	-33,3	-352,77
Candidate siRNA E2 (BTK)	43	0	-34,5	-268,37
Candidate siRNA O1 (BTK)	29	0	-30,1	-351,62
Candidate siRNA O8 (BTK)	38	-0,4	-32,5	-454,75
Candidate siRNA O5 (BTK)	38	-3,7	-32,5	-294,53
Candidate siRNA O6 (BTK)	33	-1,1	-32	-282,27
Candidate siRNA O10 (BTK)	43	-1,9	-33,3	-267,52
Candidate siRNA E3 (BTK)	29	0	-29	-267,06
Candidate siRNA O4 (BTK)	24	-1,4	-26,7	-305,91

After integrating results from machine learning predictions and thermodynamic analysis, we selected four siRNAs for experimental validation: two targeting BCL2 and two targeting BTK (one team-designed and one previously validated sequence per gene). However, one final computational step remained before ordering sequences: structural modeling using AlphaFold. Based on combined modeling results, the most promising candidates were:

• BCL2-Targeting siRNAs
  • Candidate siRNA E5 (BCL2) – ranked 1st by thermodynamic metrics and 3rd by our machine learning models
  • Candidate siRNA D1 (BCL2) – ranked 3rd by thermodynamic metrics and 1st by our machine learning models
• BTK-Targeting siRNAs
  • Candidate siRNA E1 (BTK) – ranked 1st by thermodynamic metrics
  • Candidate siRNA D2 (BTK) – ranked 2nd by thermodynamic metrics and 3rd by our machine learning models

While thermodynamic calculations provide free energy values and docking scores, they cannot reveal the three-dimensional arrangement of molecules. By modeling the Ago2–siRNA–mRNA complex with AlphaFold ¹⁶, we obtained a structural visualization that confirmed whether the siRNA guide strand aligned correctly in the Ago2 binding pocket and whether the target mRNA was properly positioned for cleavage. This step was essential to ensure that our chosen siRNAs were not only predicted to bind strongly, but also to interact in a geometrically feasible way with the silencing machinery.

We made the following assumptions while implementing our structural modeling for siRNAs:

Based on the above figures we made the following observations:

• Overall Structure
• The siRNA is positioned to some extent within the Ago2 protein, forming a duplex with the target BCL2 mRNA.
• The guide strand appears somehow oriented in the Ago2 binding channel, suggesting a plausible RISC-like conformation.
• Secondary structural elements of Ago2 surround a part of the siRNA duplex, indicating medium encapsulation.
• Functional Interpretation
• The positioning of the siRNA within Ago2 suggests correct loading, with the guide strand exposed toward the mRNA duplex which is essential for cleavage activity.
• The extended mRNA duplex conformation indicates that this design provides good complementarity and accessibility to Ago2.
• This structural outcome supports the designed BCL2 siRNA as a promising candidate for efficient silencing.
• Confidence Scores
• pTM = 0,81: suggests that the backbone arrangement of the complex is robust.
• ipTM = 0,44: AlphaFold predicts that the duplex does interact with Ago2, but its exact positioning may be dynamic or variable.

Based on the above figures we made the following observations:

• Overall Structure
• The siRNA–mRNA duplex sits firmly within the Ago2 binding channel.
• Compared to the designed siRNA, the experimentally validated one appears more extended and linear in conformation, spanning across Ago2’s central groove.
• The guide strand is positioned in direct contact with Ago2’s catalytic region, which is consistent with effective cleavage.
• Functional Interpretation
• The extended duplex conformation supports efficient recognition and alignment of the target mRNA.
• The high ipTM value indicates that AlphaFold predicts the Ago2–siRNA interaction with strong reliability, reinforcing the known experimental evidence of this siRNA’s effectiveness.
• This prediction suggests tighter anchoring and reduced conformational flexibility compared to our designed BCL2 siRNA.
• Confidence Scores
• pTM = 0,83: very high confidence in the global fold of the entire complex.
• ipTM = 0,70: strong confidence in the interface arrangement between Ago2 and the siRNA–mRNA duplex. This is significantly higher than for the designed siRNA (0.44), suggesting more reliable positioning of the duplex inside Ago2.

Based on the above figures we made the following observations:

• Overall Structure
• The BTK-targeting siRNA–mRNA duplex binds along Ago2.
• The duplex threads through the Ago2 groove but appears less tightly anchored, with one RNA end projecting more freely.
• Ago2 secondary structures engage with the duplex core.
• Functional Interpretation
• The moderate interface score (ipTM 0.39) indicates that AlphaFold predicts Ago2 binding, but is not confident in the exact positioning of the BTK siRNA–mRNA duplex. This could reflect either true biological flexibility of this siRNA within Ago2, or lower structural stability compared to the validated BCL2 siRNA.
• Despite this, the duplex still aligns broadly in the canonical RISC configuration, meaning the design is structurally plausible.
• Confidence Scores
• pTM = 0,82: high confidence in the global folding of the Ago2–RNA complex.
• ipTM = 0,39: low-to-moderate confidence in the interface positioning of the siRNA–mRNA duplex within Ago2. This suggests considerable flexibility and uncertainty in docking.

Based on the above figures we made the following observations:

• Overall Structure
• The BTK siRNA–mRNA duplex is to some extent accommodated within Ago2.
• The duplex fits more deeply in the Ago2 central channel, showing tighter contacts with structural domains.
• Functional Interpretation
• The interface confidence suggests a realistic but flexible binding mode, capturing the fact that siRNA–mRNA duplexes can shift conformation inside Ago2.
• The modeled positioning is consistent with canonical RISC loading, suggesting functional compatibility.
• Confidence Scores
• pTM = 0,82: strong confidence in the global fold, confirming the overall Ago2–RNA structure is robust.
• ipTM = 0,52: intermediate confidence in the interface alignment. This value is notably higher than that of the designed BTK siRNA (0,39)

Structural modeling of siRNA–mRNA duplexes in complex with Ago2 provided valuable insights into the comparative behavior of designed and experimentally validated candidates targeting BCL2 and BTK.

BCL2 Candidates

For BCL2, the designed candidate (D1) was predicted to adopt a plausible RISC-like configuration, with the guide strand oriented inside the Ago2 binding channel. However, its moderate interface score (ipTM = 0.44) indicates uncertainty in the precise docking of the duplex. By contrast, the experimentally validated candidate (E5) showed a much stronger interaction with Ago2 (ipTM = 0.70), with the duplex firmly anchored in the central groove and the guide strand positioned near the catalytic site. This is consistent with its known experimental efficacy, and suggests tighter anchoring and reduced conformational flexibility compared to the designed candidate.

BTK Candidates

For BTK, the designed siRNA (D2) achieved only modest interface confidence (ipTM = 0.39), reflecting either biological flexibility or limited structural stability within the Ago2 channel. In contrast, the experimentally validated candidate (E1) displayed an improved interface score (ipTM = 0.52), with the duplex positioned more deeply in Ago2’s central channel, indicating a more realistic and stable binding mode.

Global Structural Reliability

Across all models, the global structural folds (pTM ~0.82–0.83) were highly reliable, meaning the Ago2–RNA complexes themselves are robustly predicted. The main variability lies in the siRNA–Ago2 interface (ipTM values), which effectively captures the relative stability of different candidates.

Final Comments

Taken together, these results demonstrate that AlphaFold modeling was able to distinguish between designed and experimentally validated siRNAs by their predicted structural stability and docking within Ago2. The experimentally validated siRNAs (BCL2-E5 and BTK-E1) consistently showed higher interface confidence and more canonical RISC-like conformations, supporting their superior functional potential. Meanwhile, the designed siRNAs (BCL2-D1 and BTK-D2) exhibited plausible but more flexible binding modes, suggesting they may be functional but with lower predicted stability.

The targeted lipid nanoparticles play an essential role in siREN, since they are utilized for specific siRNA delivery to the CLL cells. By employing ROR1-targeting lipid nanoparticles, we take advantage of the high expression of the protein on the surface of CLL cells and the protein absence on healthy tissues, and ensure a highly-specific delivery system. Although our experiments focus on the siRNA pillar of siREN, in this page we have modeled the second pillar of our project. The modeling process gave us a better understanding of the lipid composition and membrane and set the foundations and prepared us for implementing the second pillar of siREN.

The first step in the Lipid Nanoparticle Modeling pipeline was to model and analyze all the candidate antibodies we found through literature review and select the most promising one, based on the modeling results. Our therapeutic strategy of designing a targeted LNP delivery specifically to CLL cells relies on the recognition of the ROR1 surface receptor. As we have mentioned on our Design Page, we identified three published ROR1-targeting antibodies, retrieved their Complementarity-Determining Regions (CDRs) and Framework Regions (FRs), and reconstructed their Fv fragments (VH and VL domains). More specifically, these three antibodies were ERR1-TOP43, XBR1-402 and ERR1-306 ^{24 25}. These Fv sequences served as the starting point for our computational antibody design pipeline.

The purpose of our modeling was to compare candidate antibodies based on their predicted structural compatibility and binding strength to ROR1, so that we could select the most promising binder for Fab' fragment design. This modeling process allowed us to evaluate which antibody is most suitable and promising for conjugation on the LNP surface, while avoiding the need for extensive experimental screening at this stage.

To assess binding at the molecular level, we used the HDOCK server to dock each reconstructed Fv fragment onto the extracellular domain of ROR1. This provided initial docking poses and complex structures, which were then refined for further evaluation. From HDOCK, we obtained docking scores reflecting binding energy, complementarity, and stability of the interaction.

Following docking, we used Rosetta InterfaceAnalyzer to evaluate the antibody–antigen complexes. This tool allowed us to calculate key biophysical parameters that describe the quality of the interface between the Fv fragment and the ROR1 extracellular domain. Because raw docking outputs may contain steric clashes or unrealistic conformations, we performed relaxation of the antibody-antigen complexes using Rosetta to optimize side-chain orientations and minimize structural strain.

We then used Rosetta InterfaceAnalyzer to calculate a set of quantitative parameters that describe antibody–antigen interactions, including:

• ΔG_separated (binding free energy): Estimate of binding affinity between antibody and antigen.
• dSASA_int (interface buried surface area): Indicates how much surface is buried at the interaction site, reflecting contact stability.
• hbonds_int (interface hydrogen bonds): Number of stabilizing hydrogen bonds across the interface.
• sc_value (shape complementarity): Measures geometric fit between antibody and antigen surfaces.
• nres_int (interface residues): Number of residues involved in the interaction.

These parameters allowed us to rank the three candidate Fv fragments, focusing both on energetic stability and structural complementarity. The table below includes the results we retrieved from both HDOCK Server and Rosetta InterfaceAnalyzer. All numbers were rounded up to the second decimal digit.

Parameter	Explanation	Goal	ERR1-TOP43	XBR1-402	ERR1-306
Docking Score with Ago2 (kcal/mol)	Indicates predicted binding affinity (HDOCK)	More negative is better	-245,06	-234,17	-278,74
dG_separated (REU)	Indicates the antibody-antigen binding affinity	More negative is better	-35,25	-67,56	-45,16
dSASA_int (Ų)	Indicates how much surface is buried in the interface (reflects stable contacts)	Larger values are better	1942,01	3518,15	2325,27
hbonds_int (number)	Indicates the stability of the complex	Higher is better	10	18	10
sc_value (number)	It is the shape complementarity index, values ranging from 0 to 1	Higher is better	0,53	0,58	0,61
nres_int (number)	It is the number of interface residues	Higher is better	86	150	95

Based on the computational analysis, XBR1-402 emerged as the most promising candidate for conjugation onto our lipid nanocarrier. While all three antibodies demonstrated favorable docking poses with ROR1, the quantitative evaluation highlighted XBR1-402 as the superior binder in several key aspects:

• Docking Performance
  • It achieved a favorable docking score with HDOCK (–234,17 kcal/mol), confirming a stable initial binding pose.
• Energetic Stability
  • It achieved the most favorable binding free energy (ΔG_separated = –67,56 REU), significantly lower than the other candidates, indicating stronger predicted affinity.
• Interface Contacts
  • It displayed the largest buried surface area at the interface (dSASA_int = 3518,15 Ų), reflecting more extensive molecular contacts and improved stability.
  • It formed the highest number of hydrogen bonds (18) across the interface, further stabilizing the antibody–antigen complex.
  • Its interface residue count (150) was also the largest, suggesting a broader interaction footprint.
• Comparative Evaluation
  • Although ERR1-306 showed the best shape complementarity (0,61), XBR1-402 presented a well-balanced profile across all parameters, outperforming the others in energetic and stabilizing metrics.

Taken together, these results indicate that XBR1-402 combines strong binding affinity, extensive contact formation, and structural stability, making it the most reliable choice for Fab fragment design and subsequent LNP surface conjugation.

After an extensive literature review on lipid compositions suitable for antibody conjugation, and following the guidance of Dr. Pippa, Dr. Karatasos and Dr. Demetzos, we finalized the components for each of our two proposed lipid compositions: DLin-MC3-DMA, DSPC, Cholesterol, DMG-PEG, DSPE-PEG-MAL, 18PA, DSPE-PEG.

Before advancing to the next steps of our modeling pipeline, we aimed to understand more thoroughly the molecular structure of each component by constructing their 3D representations. After retrieving 2D structural formulas from publicly available sources in .sdf or .mol2 formats, we used OpenBabel to convert them into .pdb files, the standard 3D format commonly used in the Protein Data Bank. While this conversion provided us with usable 3D structures, we aimed for greater precision and realism. Therefore, we imported the .pdb files into Avogadro2, where we applied its geometry optimization tools. This process refined bond lengths and atomic arrangements, yielding more stable and realistic molecular conformations.

Interactive Model

To maximize the insights obtained from our model, we employed CHARMM-GUI’s Membrane Builder ¹⁷ to generate a lipid bilayer representing the interface between the inner and outer layers of our nanocarrier. Our main goals were to get a better understanding of our lipid bilayers and lay the groundwork for future Lipid Nanoparticle modeling.

By employing the Bilayer Builder Tool ^{18 19 20 21 22 23}, we created 2 simulations: one for each lipid nanoparticle composition that we propose.

DLin-MC3-DMA	DSPC	Cholesterol	DMG-PEG	DSPE-PEG	DSPE-PEG-MAL	18PA
35 %	7 %	26,6 %	1,365 %	-	0,035 %	30 %

DLin-MC3-DMA	DSPC	Cholesterol	DMG-PEG	DSPE-PEG	DSPE-PEG-MAL	18PA
50 %	10 %	38,5 %	-	1,45 %	0,05 %	-

This system can be imported into GROMACS, where we can apply force field parameters and solvent conditions to simulate the nanocarrier’s environment with greater accuracy. Through these simulations, we will be able to analyze the interactions between the hydrophobic tails of the lipid molecules and evaluate the overall stability of the structure. By conducting this further analysis for both of our proposed lipid compositions, we can directly compare their properties and get a better understanding of their structure.

Our modeling work for siREN represents a complete and thorough exploration of all the major aspects of our project design. Based on modeling results, our team selected the most promising siRNAs to evaluate, compared ROR1-targeting antibodies and got a better understanding of the lipid composition and structure. Modeling played a central role in guiding our decisions and shaping the experimental direction of our project and it provided us with a deeper understanding of the mechanisms underlying siREN.

Οur modeling is complete, but there is still room for future enhancements and advancements. These may include employing Ordinary Differential Equations Systems (ODEs) to predict the target mRNA knockdown and extending our machine learning models with additional experimental datasets to improve predictive accuracy. By building on the foundation we have established, the modeling framework of siREN can evolve into a robust, generalizable platform for rational design of RNAi-based and LNP-based therapeutics.