Abstract
miRNAs are promising diagnostic biomarkers. These small RNA molecules are always present in the human body but become dysregulated when a person develops certain diseases. Although the detection of these biomarkers in cell-free tests is ongoing work, current systems often focus solely on detecting the presence or absence of a specific miRNA, rather than the miRNAs concentration. Thus, these tests may miss relative changes in miRNA concentration when disease-induced dysregulation occurs. This work, part of the WUR iGEM 2024 project (miRADAR), aimed to address this gap by incorporating an miRNA concentration-dependent threshold mechanism in a cell-free diagnostic test. In this system, continuous miRNA input concentrations need to be converted into a binary output signal, classifying the miRNA concentration as healthy (no output signal) or indicative of disease (strong output signal). To aid the experimental engineering of the test, here we use mathematical models to evaluate and assess different candidate networks. We apply a previously published multi-objective optimisation strategy to obtain designs that satisfy relevant constraints, such as low basal expression, high readout levels, and steep switching behaviour between low and high input miRNA concentrations. Models for three different biological mechanisms were compared based on their ability to generate the desired binary output signal. One approach used three-node protein networks (such as feed-forward loops), while the other two utilised RNA-based toehold systems. Overall, the toehold-mediated strand displacement systems demonstrated the most potential for experimental implementation. These systems are believed to be less burdensome in a cell-free environment, can be more readily engineered for new miRNA sequences, and showed high detection accuracy. Based on our results, we discuss how the inclusion of sequence-specific parameters could expand the design space of our mathematical models and how careful engineering of optimisation criteria is required to evaluate designs. Ultimately, our model-based study highlights that toehold-mediated strand displacement networks have the potential to be efficient miRNA detection systems for biosensing tools in the future.
1 Introduction
MicroRNAs (miRNAs) are a promising diagnostic marker and therapeutic agent. Research has identified numerous miRNAs with potential clinical applications in the detection and monitoring of neurodegenerative diseases such as Alzheimer’s and Parkinson’s disease (; ). The inhibition or activation of miRNAs involved in such diseases has been extensively studied as potential therapies (; ). These small, single-stranded molecules are present in the body to regulate transcriptional gene expression (). When carrying a disease, a person’s gene expression is differentially regulated and these changes, when compared to healthy controls, correlate with differential miRNA concentrations. This information can be utilised in a diagnostic tool, as measuring the change in concentration of disease-specific miRNAs can indicate the presence of that disease (; Wang and Zhang, 2020).
This principle formed the basis of the WUR 2024 iGEM project, miRADAR, where miRNAs were used in a cell-free diagnostic tool to help detect the neurodegenerative disease multiple sclerosis (MS) (). In this disease, immune cells attack the myelin sheaths of the nerves, leading to a range of symptoms including muscle weakness and loss of vision (). The current diagnostic procedure, involving MRI scans and lumbar punctures, is invasive and believed to provide inconclusive results for 10%–30% of the patients [Christa Benit MD, personal communication; (; )]. At present, MS has no cure, but treatments delaying the degradation of the myelin sheaths and reducing symptom progression exist (). The earlier treatment is started, the better the quality of life of the patient can be conserved, demanding a timely diagnosis of MS (; Ziemssen et al., 2016). Previous studies have found multiple miRNAs dysregulated by MS, of which hsa-miR-484 and hsa-miR-145 are examples (Figure 1, step I) (; ). This emphasises that novel miRNA-based tests could be a valuable addition to the diagnostic procedure of MS.
FIGURE 1
In the miRADAR project, the WUR 2024 iGEM team envisioned creating a cell-free blood test to aid in diagnosing complex cases of MS (). To achieve this, a simple genetic detection circuit would need to be freeze-dried on paper discs. Upon the application of an MS-positive blood sample, the presence of several MS-specific miRNAs would trigger the genetic circuit and produce a colour marker that the patient and medical professionals can then observe. If the blood sample leads to a change in colour of the system’s output, then this suggests the presence of MS-specific miRNAs.
The key of the test lies in the concentration level of the MS-specific miRNAs; the miRNAs will always be present, but their concentration can be up- or downregulated in patients (; ; ). Current developments in other miRNA-based cell-free tests either do not take this into account and focus solely on detecting the presence of a specific miRNA, or depend on a difference in visual output, which is not sensitive enough when multiple miRNAs are detected in a single test (; Wang et al., 2023). Therefore, the addition of a concentration-dependent module is an important next step. Ideally, this module would give a binary output, where either the miRNA concentration is classified as healthy or indicative of the disease (Figure 1, step II). This signal conversion could be achieved by implementing a threshold mechanism, which distinguishes whether a miRNA is below or above a threshold concentration associated with healthy patients. No system output is generated when the input miRNA concentration is considered healthy (below threshold), while a large increase in system output is generated when the input miRNA concentration is considered indicative of MS (above threshold). The sharper this switch is, the more accurate the miRNA-based diagnostic test will be. Afterwards, individual binary signals can be integrated into a modular detection module that produces a single output signal allowing for the detection of multiple different miRNAs in a single test (Figure 1, step III). Two biological mechanisms that have the potential to create the desired threshold in the dose-response curve: i) a protein-based feed-forward loop (FFL) and ii) RNA-based toehold-mediated strand displacement (TMSD) systems.
1.1 Protein-based networks: feed-forward loops
Many genetic circuits found in nature contain core network motifs consisting of a limited number of components (). One example of such a network motif are feed-forward loops (FFLs) where three nodes can interact (in)directly with each other through activation and inhibition [Figure 2A; ()]. FFL nodes can consist of an interplay between transcription factors, proteins, DNA, and RNA, and are consistently found to regulate processes like adaptation, noise filtering, biphasic behaviour, and oscillations (; ; ; Zhang et al., 2016). The basic structure of an FFL consists of a direct interaction between input node A and output node C, combined with an indirect reaction through node B. The FFL is called coherent if the direct effect of node A on node C and the indirect effect through node B are consistent. An example of a coherent network is that A activates C directly, but also activates B which also activates C in turn: in this case, A has a net positive effect on C via both paths. If the effect through both pathways is antagonistic, the FFL is considered incoherent. Next to these main reactions from node A to node C, additional regulation between the nodes is possible, including autoregulation and feedback which allows for more complex behaviour to emerge (; ).
FIGURE 2
Previous research has shown that such networks may be relevant for our miRNA detection system. While looking for networks that show adaptation using three-node networks, (
1.2 RNA-based networks: toehold-mediated strand displacement
The goal of the miRADAR project was to develop a cell-free paper-based test, which increases accessibility for patients. To create an efficient cell-free system, we are required to limit the number of biological components needed to produce an output. As we envision that our FFL systems consist of transcription and translation of node components, our cell-free test will also require compounds to process DNA and RNA into protein which will accelerate energy usage and could limit output responses (
The TMSD system has previously been integrated into miRNA diagnostic tests as an amplification module (
To improve on the initial TMSD system, (
FIGURE 3

(A) Multi-objective optimisation of the toehold-mediated strand displacement with fuel (TMSD-F) system. I) The threshold reaction: input strand anneals to the threshold strand to produce a waste strand . II) The fuel reaction: the fuel strand releases input from to produce and . III) The TMSD reaction: input strand anneals to to form and intermediate output . IV) The reporter reaction: intermediate output displaces the quencher resulting in fluorescence and a waste strand . With a longer toehold, reaction I can proceed at a faster kinetic rate than reactions II-IV (all ). The universal toehold is shown in pink. Image adapted from (
The miRADAR project also considered a further simplification of the TMSD-F system (
FIGURE 4

(A) The TMSD-NF system analysed by the WUR 2024 iGEM team. I) The threshold reaction: input anneals to threshold strand to form a complex . II) The TMSD reaction: binds to to release intermediate output and . III) Aptamer reaction: aptamer anneals to intermediate output to produce fluorescence (denoted by a star). The toehold is shown in orange. (B) Normalised dose-response curves for the FFL (blue), TMSD-F (orange), TMSD-NF systems (pink) and a perfect binary response (purple). The curves were normalised by dividing the fluorescent output for each dose by the maximum system output. = 50 nM. (C) Parameter space exploration of the TMSD-NF system. The ratio (log scale) is plotted against the score (maximum of 10). The colour bar represents the .
1.3 Multi-objective optimisation algorithms for model design
FFL and TMSD are biological mechanisms that could, with further optimisation, form a sharp threshold in the dose-response curve. With the described kinetic models, we wish to optimise the parameters of our systems to increase the accuracy and steepness of the threshold while additionally lowering the basal expression before that threshold point (Figure 1, step II, see a schematic of these behaviours labelled -). Previous research has varied parameters of FFL networks to observe how parameter relationships impact FFL performance (
In multi-objective optimisation problems, trade-offs between different objectives are likely. As an example, one could envision that, for the threshold mechanisms, a higher slope (i.e., the sharpness of the increase between “healthy” and “disease” in Figure 1, step II, ) would reduce the threshold accuracy (i.e., is the lowest output-producing input concentration close to the desired switching/threshold concentration; in Figure 1, step II). The algorithm from (
By applying the above-described optimisation strategy of (
2 Methods
2.1 Mathematical models
2.1.1 Three-node networks and feed forward loops
The feed forward loop (FFL) system consists of three nodes, which can inhibit or activate each other and themselves (Figure 2A). We assume these connections to represent transcription and translation where one node’s mRNA is translated to a transcription factor that can regulate preceeding nodes. A single node connection, denoted , in our three-node networks can be described by combining the regulation type and regulation strength into , where node acts on node (
The ODEs were solved with a CVODES solver provided by Serban and Hindmash which was adjusted by
2.1.2 Fuel-regulated toehold mediated strand displacement system
The ODEs for toehold mediated strand displacement system with a fuel reaction (TMSD-F) were derived by
The ODEs were solved with the MATLAB CVODES stiff solver under standard parameters, except absolute tolerance = and relative tolerance = so solutions with high accuracy could be achieved (
2.1.3 Toehold mediated strand displacement systems without fuel
The ODEs for the fuel-removed TMSD system (TMSD-NF) were adapted from TMSD-F and consist of an equation for the threshold reaction (with rate ), the TMSD reaction (with rate ) and the binding to a fluorescent aptamer (with rate ; Figure 4A). The complete set of model reactions and equations are in Supplementary Methods S1.5. The ODEs describe the change over time of the input miRNA , the threshold sequence , the resulting complexes formed with the intermediate output and the output sequence .
The ODEs were solved with the MATLAB CVODES non-stiff solver using standard parameters, except absolute tolerance = and relative tolerance = (
2.2 Mathematical representation of a dose-response curve
Mathematically, dose-response curves can be described with the Hill functionwhere is the concentration of input miRNA that results in half of the differential system output (
Simulating a large number of doses would slow down the optimisation drastically, so a minimal amount of doses, which still accurately represent the curve, had to be determined (Supplementary Methods S1.1). In total, 19 doses are necessary to create an accurate dose-response curve, with 10 of the total doses centred around the threshold dose. This is shown schematically in Figure 1 (step II).
2.3 Optimisation objectives
The objectives of the optimisation of the FFL and TMSD-F systems are based on the desired dose-response curve (Figure 1, step II). Ideally, this curve has i) a high threshold accuracy, ii) a high system output after the threshold is passed, iii) a low basal expression before the threshold, and iv) a steep slope.
The from the Hill function (Equation 1) should be at the dose where the miRNA-dependent threshold should be. Therefore we want to minimise the difference between the threshold concentration we want with the threshold concentration we obtain with our model . We define our first objective to minimise this difference followingwhere is calculated with the Hill function obtained from an optimisation solution and is the desired threshold value. This objective was minimised in the optimisation, and we will refer to this as minimising the threshold differences or maximising threshold accuracy. The other three objectives were set as optimisation constraints, meaning that objectives , and are constrained to user-input values that comply with well-performing dose-response curves.
The second objective, which assesses the maximum system output can be defined as constraintwhere the lower boundary and the upper boundary were set to reach the known maximum system output. For FFL the boundaries were set to and , while for TMSD-F . By defining then we are enforcing our output at , , to be greater or equal to the output produced at the lowest input concentration, . Next, the basal expression was quantified by summing up the system outputs induced by the four lowest doses of input miRNAwhere represents the index of the dose and both and are a user-input value to constrain the basal expression. In this work, is always set to 0, while is set to 2 for the FFL optimisation and 7.5 or 15 for the TMSD-F optimisation. The decision to use the four lowest input doses is discussed in Supplementary Methods S1.1, 1.2 where we show that simulating more input doses increases computational time to optimise a single model, whilst including more doses within this scoring function does not influence the score per data point for an optimal model.
According to (
In the optimisation of FFL and three-node networks, a fifth constraint objective was added to limit the number of node connections as followswith and resulting in a summation of 9 absolute (regulation effect) values in total. The bounds and are user-input integers and were defined here as and to limit the total connections in our three-node system to a maximum of 4. In a three-node network the maximum amount of node connections is 9 and we do not enforce the FFL structure on our networks.
2.4 Multi-objective optimisation algorithm
The optimisation algorithm combined the global solver eSS (enhanced scatter search; release 2010A) and the local solver misqp (for FFLs and three-node networks; version 7.1) or fmincon (for TMSD systems; version MATLAB 2024a) into a hybrid solver (
2.5 Latin hypercube sampling for parameter space exploration of TMSD-NF
The TMSD-NF system consists of three reactions rates (, and ). This allows us to more extensively assess the relationship between these parameters and system performance. To do this, we sampled permissible parameter space for the TMSD-NF system with a Latin Hypercube sample ( = 10,000). The objectives presented above were slightly adjusted to better fit the behaviour of the TMSD-NF system, which has linear dose-response curves that cannot be described by the Hill function.
The maximum output of the system is determined by the known concentration of aptamer (Supplementary Methods S2.3.1). Therefore, we know what the expected maximum output of the system should be . We, therefore, scored how close the observed maximum output is to the expected maximum asto score the basal expression before the threshold, the formula remained unadjusted from above, except the lower and upper boundaries were removed:to measure the steepness of TMSD-NF systems, we replaced function (Equation 5) withthis change was required since TMSD-NF systems produced optimal dose-response curves with sharp transitions between healthy and disease regimes (Figure 1, step II) once the threshold input dose has been crossed - a qualitatively different behaviour to which we had before. This means that fitting a Hill function, approximating and using as a scoring metric, in these cases became unreliable in high-throughput since the concentration where the output was half the maximal value could not be fixed to . In the Supplementary Methods S1.2 we discuss the impact of this change and show that, for simulated test cases, the formulation of produces qualitatively the same results as function .
The change in gradient just before and after indicates the threshold accuracy, where a higher value correlates with minimal differences between the observed threshold and the expected threshold. The corresponding scoring function is defined as
The final score was computed aswhere and and . The lower boundary of the rescaling was set to 0 and the upper boundary of the rescaling was set to 10. and are the minimum and maximum values for each corresponding set of scoring function outcomes from our sampling. The rescaling to a standardised range of 0–10 was necessary for a fair comparability of the different scoring functions as they originally varied over different ranges. Without this, scoring functions with large outcome values would disproportionately influence the final score. The scoring function was reverse-coded to convert the lowest basal expressions to the highest scores. The scoring function was assigned a lower weight because it was consistently observed to achieve satisfactory values and thus less helpful in differentiating the curves.
3 Results
3.1 Optimising three-node networks and feed-forward loops
The first mechanism applied to convert continuous input concentrations of miRNA into a binary fluorescent signal are three-node networks, of which the FFL system is a special example (Figure 2A). In the multi-objective optimisation process (Section 2.3, Equations 2 - 6), the threshold accuracy , i.e. how close the of the dose-response curve is to , was optimised under additional constraints on the basal expression , slope and the number of node connections in the network . The closer the value for is to 0, the smaller the difference is between the system’s value of and what we wish to achieve. Conversely, this could be a considered as maximising the threshold accuracy. We visualise our results over a two-dimensional - sample space, i.e. the threshold difference is plotted against the slope from the estimated Hill function. Multiple intervals of slope values , ranging from 0 to 20 in increments of 2, were tested, resulting in ten combinations of threshold accuracies and slope values . Based on our results, we observed a trade-off where higher slope values combine with lower threshold differences . The effect of limiting the basal expression or the number of node connections in our three-node networks was judged according to their effect on and .
According to the plotted search space, our three-node systems consistently reached minimal threshold differences, showing that high threshold accuracies are robust to changes in the slope of output dose-response Hill functions (Figure 2B; Supplementary Figure S4). These values were reached regardless of the limits set by and , indicating that low basal expressions and fewer node connections do not compromise the slope and threshold accuracy. In other words, it is possible to construct a three-node system complying with all the set objectives. Outliers with poor threshold accuracies exist but this was an issue for every constraint combination, suggesting that the algorithm might sometimes be stuck in a local minima (Supplementary Figure S4).
Based solely on the objective values, multiple well-performing systems exist but this is not reflected in the simulated dose-response curves. There, some curves show irregular behaviour, where the output at higher doses of input miRNA is not constant. As an example, from two systems performing similarly in our - search space, one dose-response curve shows irregular behaviour (Figures 2B,C, dagger), while the other does not (Figures 2B,C, asterisk). Their respective responses over time expose that systems with irregular dose-response curves do not reach steady states, but instead, oscillate (Figure 2D). The differences in behaviour are reflected in the topologies of the system (Supplementary Figure S6). Both non-oscillating and oscillating FFL mechanisms involve strong negative regulation on A (either through , or ). Our non-oscillating systems tend to form more regular FFL systems or linear pathways, but oscillating systems form a Goodwin oscillator [positive , positive and negative ; (
3.2 Optimising fuel-regulated toehold systems
The TMSD system solely relies on nucleotide binding, providing an advantage over the energy-demanding transcription and translation necessary for the functioning of an FFL. The TMSD-F system contains a fuel strand, which catalytically speeds up and increases the production of fluorescent output (Figure 3A II). The same objectives as for the three-node networks, except for , were applied to optimise the reaction rates and improve the dose-response curve produced with the TMSD-F system from (
Basal expression limits (, Equation 4) of 7.5 nM and 15 nM, alongside unconstrained , were used as constraints for TMSD-F optimisation. The effect on the threshold accuracy (, Equation 2) and slope (, Equation 5) were again evaluated over a two-dimensional space. A restriction of to 15 nM did not affect the threshold accuracy of the switch compared to unconstrained basal expression (overlap of orange and blue dots in Figure 3B). However, further reducing the basal expression to 7.5 nM detrimentally reduced the threshold accuracy, especially for lower slope values (red dots in Figure 3B). The plotted dose-response curves revealed that the curve shifts slightly to the right of the when was limited to 7.5 nM, which reduced the threshold accuracy (Supplementary Figure S7). In contrast to our three-node networks and FFL systems, no oscillations are observed in the dose-response curve.
From the observed search space, we discovered that constraining our slope to be between 16 and 18 and the basal expression to be at most 15 nM provided a desirable solution for our purposes. The resulting network combines the minimal difference between simulated and expected value, a high slope , and a limit on the basal expression (Figure 3B, asterisk). Compared to (
With these results in mind, the improved threshold accuracy, higher slope and reduced basal expression of the optimised TMSD-F produces the dose-response curve shown in Figure 3C. Compared to our three-node systems, the better engineering possibilities and the absence of oscillations are great advantages for employing the TMSD-F system as the concentration-dependent module in our miRNA diagnostic test.
3.3 Optimising toehold systems in the absence of fuel reactions
As the kinetic model of the TMSD-NF system was adapted from TMSD-F, we assumed that the rates of similar reactions could be transferred between the systems (Figure 4A). Therefore, the threshold reaction proceeds again with rate , while the TMSD and aptamer binding reaction proceed with rates (here named and , respectively).
In simulations with these reaction rates, the TMSD-NF system generates more linear dose-response curves than the TMSD-F system. The TMSD-NF system showed minimal basal expression, but the slope of the dose-response curve was less steep than in the TMSD-F system (Figure 4B, pink line). To test whether the decrease in basal expression was due to removing the fuel reaction, the initial concentration of the fuel component was set to 0 in the TMSD-F system, thereby eliminating the fuel reaction (Supplementary Figure S12). The simulations showed a small increase of basal expression in the absence of the fuel reaction.
To further engineer the TMSD-NF system, the initial concentrations of intermediate output (in complex with ) and aptamer were adjusted to observe the trade-off between maximum system output and the slope of the resulting dose-response curve (Supplementary Figures S17, S18). Concentrations of and equal to produced dose response curves with steep slopes and high levels of fluorescence when miRNA inputs are above (set to 2 nM in Supplementary Figure S18). Consequently, though, the fluorescent output cannot exceed , which could be problematic when trying to detect miRNA with low threshold concentrations between off- and on-states.
In
To test this observation, a Latin Hypercube sample (Section 2.5) was created of the permissible parameter space to get a better understanding of the relationship between these TMSD-NF parameters. The permissible parameter space is defined by scoring functions s1 to s4 (Equations 7 - 10). Upon plotting the ratio against the performance score (Equation 11), a clear trend emerged. The larger the ratio, the better the dose-response curve is observed by high score values. From roughly onwards, the TMSD-NF system always manages to produce the intended outcome (Figure 4C). Interestingly, has a defining character when the ratio is small. Even though , a low-binding aptamer (low ) led to a high score. Potentially, the output TMSD strand is produced much quicker, but must accumulate to larger quantities to overcome the poor aptamer binding imposed by . In this sense, the aptamer creates the threshold instead of the antisense strand . Nonetheless, our results suggest that a threshold reaction should have a that is at least 1,000 times faster than . This can be achieved by extending the toehold length of the threshold reaction, thereby increasing , or shortening the toehold length of the TMSD reaction described by .
3.4 Comparison of the three systems
As many miRNAs can differentiate between healthy people and patients diagnosed with MS, we want to create a system that can respond to a suite of differentially-expressed miRNAs. However, each miRNA will have a different threshold concentration that distinguishes between healthy patients and those with MS. Therefore, the concentration-dependent module requires a modular design that is easily adaptable to new values.
By changing the input concentrations of the system, where the threshold strand concentration is equal to , we can study how adaptive the systems are to new threshold values (Figure 5). The switching behaviour of the TMSD systems (yellow and pink lines) track the value of used in the simulations, indicating good adaptive behaviour to new input miRNAs necessary for a modular system. Conversely, regardless of input miRNA’s values, the FFL system will always produce the same dose-response curve, which converts the miRNA input to a binary output at equal to 50. For a new threshold value, the optimisation of our three-node networks would have to be redone to find the best topology, which is disadvantageous relative to the easy changes in input concentrations that can be made for the TMSD systems.
FIGURE 5

Normalised dose-response curves for FFL (blue), TMSD-F (orange), TMSD-NF (pink) and a perferct binary response (purple) at multiple values of . The curves were normalised by dividing the fluorescent output for each dose by the maximum system output. (A) = 1 nM. (B) = 25 nM. (C) = 50 nM. (D) = 75 nM.
Focussing back on a of 50 nM, the FFL system does have a very high threshold accuracy, an acceptable basal expression, and a sufficiently steep slope (Figure 5C). The FFL topology does not oscillate over time resulting in the maximum output of the system, 20 nM, being consistent after the mark (Supplementary Figure S5). In comparison, the TMSD-NF system has an even lower basal expression and a similar slope, but the maximum output of the system is limited by the value as discussed above. The incorporation of the fuel reaction in TMSD-F introduces both benefits and limitations compared to the other two systems. The maximum possible fluorescence is the highest of the three systems and can be engineered to be even higher, but this comes at the cost of high basal expression (Figure 5C). Limiting the basal expression further might further slow down the system which would be undesirable depending on the system’s application. The TMSD-F system does accurately detect the threshold, and the slope in the dose-response curve is the steepest out of the three systems.
In the final diagnostic test, our results recommend the use of TMSD systems as the chance of oscillations and the limited scalability of three-node networks, like the FFL system, are undesired. The TMSD-F system has the most potential if further reduction of the basal expression can be achieved. TMSD-F has a high threshold accuracy, the highest slope and the highest maximum system output. The latter point, in particular, is a disadvantage of the TMSD-NF system, where the maximum system output is relatively low. If further optimisation of the TMSD-F system proves difficult, the system would be best employed at low values as the increase in system output provided by the fuel is most relevant at these input concentrations. At higher values, the basal expression of the TMSD-F system can increase and cause false positives, so for those values, the TMSD-NF system might be preferred. Ultimately, we have shown the TMSD systems have the theoretical potential to be simply engineered for application as concentration-dependent modules in a range of miRNA-based detection tools.
4 Discussion
In this work, we have utilised a previously published multi-objective optimisation strategy to design biological mechanisms that are capable of converting (continuous) miRNA inputs into binary output signals. As per the last section of the results, the RNA-based TMSD systems outperform the protein-based three node FFL system. These TMSD systems can easily be adapted to other input miRNAs (with different values), have a high threshold accuracy, and do not cause aberrant behaviour such as oscillations. Crucially, though, our key insights into the functioning of toehold systems are yet to be experimentally validated. The previous success whereby models of toehold systems have been experimentally validated (see the Supplementary Material of
4.1 TMSD engineering
The fuel reaction in particular is an interesting target for optimising the TMSD-F system further. Our results showed that this reaction in the TMSD-F system results in a trade-off between a high fold-change in the system output versus lower basal expression (Figure 3B). The optimisation results showed that the basal expression can be decreased by lowering the rate describing the fuel, TMSD, and reporter reaction (Figure 3A, reactions II, III and IV, respectively). Biologically, this could be achieved by shortening the toehold length which initiates the reaction (Zhang and Winfree, 2009;
A critical limitation of the current TMSD-F and TMSD-NF models is the systems’ reliance on domain binding (i.e., toehold to toehold) rather than sequence-specific binding. While, for example, the length of the toehold is essential for the speed of the TMSD reaction rates, the sequence itself also plays a role (Zhang and Winfree, 2009;
Including sequence-specificity in the model becomes even more evident when considering the application of the TMSD system in diagnostic tests. Ultimately, our designed system would be used in applications to detect multiple miRNAs simultaneously, meaning that multiple TMSD systems will need to work in parallel. Here, sequence specificity becomes crucial, as the wrong miRNA should not trigger a TMSD reaction and produce false positives or negatives. In this work, we assumed that parallel detection is possible, allowing us to model one TMSD system that can be repurposed for all miRNA that we wish to detect. This assumption could potentially be violated on sequence level, which could, for example, lead to a decoy miRNA with a slight mismatch to bind to the threshold strand of the target miRNA. This could cause false positives where the concentration of the target miRNA did not pass the threshold but, together with the decoy miRNA, the threshold is surpassed. This signifies the need for well-designed toeholds that are highly specific for one miRNA only. Fortunately, when TMSD was used as an amplification module, it was specific to single nucleotide mismatches (Zhang et al., 2020). Other work underlines the importance of sequence specificity, but current models of this mismatch effect are dependent on specific toehold lengths and the position of the mismatching nucleotide (
To tackle these issues of sequence specificity and decoy miRNAs for our designed systems, we propose three extensions to our work for practitioners and future research through the use of Figure 4C. In the first instance, our modelling framework could be extended to incorporate sequence specificity by making use of the previously developed KinDa tool. This tool compares the functioning of TMSD systems at the domain and sequence level with stochastic modelling (
Alternatively, extra experimental data could be obtained to further evaluate the performance of TMSD systems. For example, in the first instance, practitioners could evaluate the performance of the TMSD system with varying amounts of initial aptamer concentration or testing aptamers of different binding strength. As we observe in Figure 4C, when the parameter (related to concentration of aptamer) is varied then at low ratios we observe varying TMSD performance. Consequently, if TMSD performance is shown to depend on aptamer concentration or sequence then this correlates with a low ratio suggesting there is an issue within the system (e.g. decoy miRNAs or nucleotide mismatches could be present in the sample perturbing the TMSD). Second, the robustness of a TMSD system’s response to a particular miRNA could be evaluated in conjunction to dose-response curves obtained when using decoy miRNA or miRNA with mismatching nucleotides as inputs. In the event that significant overlap of dose-response curves is observed between these conditions, then this is suggestive of the TMSD system being insensitive to changes in nucleotide sequence since mismatching miRNA inputs can trigger the TMSD system as well as perfectly matching miRNA inputs. In both instances — either when TMSD systems are sensitive to aptamer alterations or TMSD performance significantly overlaps between perfect and mismatching miRNA inputs — then we would encourage testing other TMSD designs for other potential biomarkers to find robust detection mechanisms.
4.2 Improving optimisation for better design of FFL threshold mechanisms
A major issue in our three-node network designs is the formation of topologies that cause oscillations over time in the fluorescent output. For a correctly working threshold mechanism, the system should reach steady state within a reasonable time period. Otherwise, the output signal is inconsistent, and accurate measurements of the miRNA concentrations are difficult. Although the topologies of the networks causing oscillations and the networks resulting in smooth dose-response curves do not entirely overlap, they share heavy negative regulation on node (Figure 4A; Supplementary Figure S6). In previous studies, this type of negative feedback has been associated with networks that provide robustness to noise, as well as oscillations (
However, strong negative feedback is also associated with oscillating networks (
Therefore, adding a constraint to the optimisation strategy that prevents any solutions with oscillations would be necessary. This could be achieved by adding careful constraints to which reactions within a network are allowed, and is required since our current constraints are insufficient to achieve this currently. Alternatively, the method of
If we take a step back and evaluate the multi-objective optimisation framework as a whole, we have observed trade-offs between different objectives through the visualisation of our search spaces in Figures 2B, 3B. However, decisions cannot be based on this information alone. Careful examination of the dose-response curves and time-responses was necessary to filter out undesired behaviour, like oscillations produced by three-node networks, and determine the influence of a smaller on the time required to reach the steady state in the TMSD-F optimisation. With additional constraints, these filtering steps could be included in future design strategies based on the approach used here.
5 Conclusion
In summary, this study modelled and explored three biological mechanisms in their ability to convert continuous miRNA input into a binary output above a specific threshold. All the system designs studied here showed potential for future use in sensor- or diagnostic tests. However, the RNA-based TMSD systems are easier to engineer, more stable, and more adaptable to new input miRNAs than protein-based networks such as the FFL system. The TMSD-F system would outcompete the TMSD-NF system at higher threshold values if the basal expression produced by our TMSD-F design could be further reduced.
In the future, the miRADAR project of the WUR iGEM 2024 team envisions the incorporation of this concentration-dependent module into cell-free miRNA diagnostic tests (
Statements
Data availability statement
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://git.wur.nl/ssb/publications/designing-mirna-detection-networks.
Author contributions
RV: Conceptualization, Writing – original draft, Visualization, Software, Writing – review and editing, Formal Analysis. RS: Formal analysis, Visualisation, Writing – original draft, Conceptualization, Supervision, Writing – review and editing.
Funding
The author(s) declare that financial support was received for the research and/or publication of this article. RS is supported by NWO-M project 17336 (“Systems biology analysis of infection structure development in a plant pathogenic fungus”).
Acknowledgments
We acknowledge the support of the WUR 2024 iGEM team, miRADAR, and thank all team members and supervisors. In particular, we acknowledge the support of the Laboratories of Systems and Synthetic Biology, Microbiology (Nico Claassens) and Bioprocess Engineering (Mark Bisschops) at Wageningen UR in coordinating the iGEM team. We are also grateful to discussions with Bram de Jonge and Pieter Candry to help us model the TMSD-NF system. Finally, we thank Maria Suarez Diez and Edoardo Saccenti for critical reading of this manuscript before submission.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Generative AI statement
The author(s) declare that Generative AI was used in the creation of this manuscript. ChatGPT with GPT version 4o was used sparingly for debugging code and for checking grammatical and spelling errors.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fsysb.2025.1601854/full#supplementary-material
References
1
AkayA.ReddyH. N.GallowayR.KozyraJ.JacksonA. W. (2024). Predicting DNA toehold-mediated strand displacement rate constants using a DNA-BERT transformer deep learning model. Heliyon10, e28443. 10.1016/j.heliyon.2024.e28443
2
AlonU. (2007). Network motifs: theory and experimental approaches. Nat. Rev. Genet.8, 450–461. 10.1038/nrg2102
3
AngJ.HarrisE.HusseyB. J.KilR.McMillenD. R. (2013). Tuning response curves for synthetic biology. ACS Synth. Biol.2, 547–567. 10.1021/sb4000564
4
BaumK.PolitiA. Z.KofahlB.SteuerR.WolfJ. (2016). Feedback, mass conservation and reaction kinetics impact the robustness of cellular oscillations. PLOS Comput. Biol.12, e1005298. 10.1371/journal.pcbi.1005298
5
BerleantJ.BerlindC.BadeltS.DannenbergF.SchaefferJ.WinfreeE. (2018). Automated sequence-level analysis of kinetics and thermodynamics for domain-level DNA strand-displacement systems. J. R. Soc. Interface15, 20180107. 10.1098/rsif.2018.0107
6
BoadaY.Reynoso-MezaG.PicóJ.VignoniA. (2016). Multi-objective optimization framework to obtain model-based guidelines for tuning biological synthetic devices: an adaptive network case. BMC Syst. Biol.10, 27. 10.1186/s12918-016-0269-0
7
ChenC.WenJ.WenZ.SongS.ShiX. (2023). DNA strand displacement based computational systems and their applications. Front. Genet.14, 1120791. 10.3389/fgene.2023.1120791
8
DoroszkiewiczJ.GroblewskaM.MroczkoB. (2022). Molecular biomarkers and their implications for the early diagnosis of selected neurodegenerative diseases. Int. J. Mol. Sci.23, 4610. 10.3390/ijms23094610
9
DuY.XieL.LiuJ.WangY.XuY.WangS. (2014). Multi-objective optimization of reverse osmosis networks by lexicographic optimization and augmented epsilon constraint method. Desalination333, 66–81. 10.1016/j.desal.2013.10.028
10
EgeaJ. A.Balsa-CantoE.GarcíaM. S. G.BangaJ. R. (2009). Dynamic optimization of nonlinear processes with an enhanced scatter search method. Industrial and Eng. Chem. Res.48, 4388–4401. 10.1021/ie801717t
11
ElmiZ.LiB.LiangB.LauY.Borowska-StefańskaM.WiśniewskiS.et al (2023). An epsilon-constraint-based exact multi-objective optimization approach for the ship schedule recovery problem in liner shipping. Comput. and Industrial Eng.183, 109472. 10.1016/j.cie.2023.109472
12
ExlerO.SchittkowskiK. (2007). A trust region SQP algorithm for mixed-integer nonlinear programming. Optim. Lett.1, 269–280. 10.1007/s11590-006-0026-1
13
GardnerD. J.ReynoldsD. R.WoodwardC. S.BalosC. J. (2022). Enabling new flexibility in the SUNDIALS suite of nonlinear and differential/algebraic equation solvers. ACM Trans. Math. Softw.48 (31), 1–24. 10.1145/3539801
14
GentileG.MorelloG.La CognataV.GuarnacciaM.ConfortiF. L.CavallaroS. (2022). Dysregulated miRNAs as biomarkers and therapeutical targets in neurodegenerative diseases. J. Personalized Med.12, 770. 10.3390/jpm12050770
15
GhasemiN.RazaviS.NikzadE. (2017). Multiple sclerosis: pathogenesis, symptoms, diagnoses and cell-based therapy. Cell J. (Yakhteh)19, 1–10. 10.22074/cellj.2016.4867
16
GiovannoniG.ButzkuevenH.Dhib-JalbutS.HobartJ.KobeltG.PepperG.et al (2016). Brain health: time matters in multiple sclerosis. Multiple Scler. Relat. Disord.9, S5–S48. 10.1016/j.msard.2016.07.003
17
GonzeD.BernardS.WaltermannC.KramerA.HerzelH. (2005). Spontaneous synchronization of coupled circadian oscillators. Biophysical J.89, 120–129. 10.1529/biophysj.104.058388
18
GoodwinB. C. (1965). Oscillatory behavior in enzymatic control processes. Adv. Enzyme Regul.3, 425–438. 10.1016/0065-2571(65)90067-1
19
HauserS. L.CreeB. A. C. (2020). Treatment of multiple sclerosis: a review. Am. J. Med.133, 1380–1390.e2. 10.1016/j.amjmed.2020.05.049
20
HindmarshA. C.BrownP. N.GrantK. E.LeeS. L.SerbanR.ShumakerD. E.et al (2005). SUNDIALS: suite of nonlinear and differential/algebraic equation solvers. ACM Trans. Math. Softw.31, 363–396. 10.1145/1089014.1089020
21
HoP. T. B.ClarkI. M.LeL. T. T. (2022). MicroRNA-Based diagnosis and therapy. Int. J. Mol. Sci.23, 7167. 10.3390/ijms23137167
22
HolehouseJ.CaoZ.GrimaR. (2020). Stochastic modeling of autoregulatory genetic feedback loops: a review and comparative study. Biophysical J.118, 1517–1525. 10.1016/j.bpj.2020.02.016
23
JungJ. K.ArchuletaC. M.AlamK. K.LucksJ. B. (2022). Programming cell-free biosensors with DNA strand displacement circuits. Nat. Chem. Biol.18, 385–393. 10.1038/s41589-021-00962-9
24
KholodenkoB. N. (2006). Cell-signalling dynamics in time and space. Nat. Rev. Mol. Cell Biol.7, 165–176. 10.1038/nrm1838
25
KimD.KwonY. K.ChoK. H. (2008). The biphasic behavior of incoherent feed-forward loops in biomolecular regulatory networks. BioEssays30, 1204–1211. 10.1002/bies.20839
26
LeeJ.NaH. K.LeeS.KimW. K. (2021). Advanced graphene oxide-based paper sensor for colorimetric detection of miRNA. Microchim. Acta189, 35. 10.1007/s00604-021-05140-1
27
LeverM.LimH. S.KrugerP.NguyenJ.TrendelN.Abu-ShahE.et al (2016). Architecture of a minimal signaling pathway explains the T-cell response to a 1 million-fold variation in antigen affinity and dose. Proc. Natl. Acad. Sci.113, E6630–E6638. 10.1073/pnas.1608820113
28
LiY.ZhouL.NiW.LuoQ.ZhuC.WuY. (2019). Portable and field-ready detection of circulating MicroRNAs with paper-based bioluminescent sensing and isothermal amplification. Anal. Chem.91, 14838–14841. 10.1021/acs.analchem.9b04422
29
LiangK.WangH.LiP.ZhuY.LiuJ.TangB. (2020). Detection of microRNAs using toehold-initiated rolling circle amplification and fluorescence resonance energy transfer. Talanta207, 120285. 10.1016/j.talanta.2019.120285
30
LuJ.ClarkA. G. (2012). Impact of microRNA regulation on variation in human gene expression. Genome Res.22, 1243–1254. 10.1101/gr.132514.111
31
MaW.TrusinaA.El-SamadH.LimW. A.TangC. (2009). Defining network topologies that can achieve biochemical adaptation. Cell138, 760–773. 10.1016/j.cell.2009.06.013
32
MachinekR. R. F.OuldridgeT. E.HaleyN. E. C.BathJ.TurberfieldA. J. (2014). Programmable energy landscapes for kinetic control of DNA strand displacement. Nat. Commun.5, 5324. 10.1038/ncomms6324
33
ManganS.AlonU. (2003). Structure and function of the feed-forward loop network motif. Proc. Natl. Acad. Sci.100, 11980–11985. 10.1073/pnas.2133841100
34
Marquez-LagoT. T.StellingJ. (2010). Counter-intuitive stochastic behavior of simple gene circuits with negative feedback. Biophysical J.98, 1742–1750. 10.1016/j.bpj.2010.01.018
35
MavrotasG. (2009). Effective implementation of the epsilon-constraint method in multi-objective mathematical programming problems. Appl. Math. Comput.213, 455–465. 10.1016/j.amc.2009.03.037
36
MiloR.Shen-OrrS.ItzkovitzS.KashtanN.ChklovskiiD.AlonU. (2002). Network motifs: simple building blocks of complex networks. Science298, 824–827. 10.1126/science.298.5594.824
37
iGEM (2024). miRADAR: WUR 2024 iGEM. Available online at: 2024.igem.wiki/wageningenur/.
38
NguyenT. P. N.KumarM.FedeleE.BonannoG.BonifacinoT. (2022). MicroRNA alteration, application as biomarkers, and therapeutic approaches in neurodegenerative diseases. Int. J. Mol. Sci.23, 4718. 10.3390/ijms23094718
39
Otero-MurasI.BangaJ. R. (2016). Design principles of biological oscillators through optimization: forward and reverse analysis. PLOS ONE11, e0166867. 10.1371/journal.pone.0166867
40
Otero-MurasI.BangaJ. R. (2017). Automated design framework for synthetic biology exploiting pareto optimality. ACS Synth. Biol.6, 1180–1193. 10.1021/acssynbio.6b00306
41
PietersP. A.NathaliaB. L.van der LindenA. J.YinP.KimJ.HuckW. T. S.et al (2021). Cell-free characterization of coherent feed-forward loop-based synthetic genetic circuits. ACS Synth. Biol.10, 1406–1416. 10.1021/acssynbio.1c00024
42
QianL.WinfreeE. (2011). Scaling up digital circuit computation with DNA strand displacement cascades. Science332, 1196–1201. 10.1126/science.1200520
43
RahmanA.TiwariA.NarulaJ.HicklingT. (2018). Importance of feedback and feedforward loops to adaptive immune response modeling. CPT Pharmacometrics and Syst. Pharmacol.7, 621–628. 10.1002/psp4.12352
44
RegevK.HealyB. C.PaulA.Diaz-CruzC.MazzolaM. A.RahejaR.et al (2018). Identification of MS-specific serum miRNAs in an international multicenter study. Neurology Neuroimmunol. and Neuroinflammation5, e491. 10.1212/NXI.0000000000000491
45
SeeligG.SoloveichikD.ZhangD. Y.WinfreeE. (2006). Enzyme-free nucleic acid logic circuits. Science314, 1585–1588. 10.1126/science.1132493
46
SerbanR.HindmarshA. C. (2005). “CVODES: the sensitivity-enabled ODE solver in SUNDIALS,” in Volume 6: 5Th international conference on multibody systems, nonlinear dynamics, and control, parts A, B, and C. Long Beach, California, USA: ASMEDC, 257–269. 10.1115/DETC2005-85597
47
SøndergaardH. B.HesseD.KrakauerM.SørensenP. S.SellebjergF. (2013). Differential microRNA expression in blood in multiple sclerosis. Multiple Scler. J.19, 1849–1857. 10.1177/1352458513490542
48
SongM.PanK.SuH.ZhangL.MaJ.LiJ.et al (2012). Identification of serum MicroRNAs as novel non-invasive biomarkers for early detection of gastric cancer. PLoS ONE7, e33608. 10.1371/journal.pone.0033608
49
StögbauerT.WindhagerL.ZimmerR.RädlerJ. O. (2012). Experiment and mathematical modeling of gene expression dynamics in a cell-free system. Integr. Biol.4, 494–501. 10.1039/c2ib00102k
50
SzekelyP.SheftelH.MayoA.AlonU. (2013). Evolutionary tradeoffs between economy and effectiveness in biological homeostasis systems. PLOS Comput. Biol.9, e1003163. 10.1371/journal.pcbi.1003163
51
TanedaA. (2015). Multi-objective optimization for RNA design with multiple target secondary structures. BMC Bioinforma.16, 280. 10.1186/s12859-015-0706-x
52
TullmanM. J. (2013). Overview of the epidemiology, diagnosis, and disease progression associated with multiple sclerosis. Am. J. Manag. Care, 19, 15–20.
53
TysonJ. J.NovákB. (2010). Functional motifs in biochemical reaction networks. Annu. Rev. Phys. Chem.61, 219–240. 10.1146/annurev.physchem.012809.103457
54
UllnerE.BucetaJ.Díez-NogueraA.García-OjalvoJ. (2009). Noise-induced coherence in multicellular circadian clocks. Biophysical J.96, 3573–3581. 10.1016/j.bpj.2009.02.031
55
WalgraveH.ZhouL.De StrooperB.SaltaE. (2021). The promise of microRNA-based therapies in Alzheimer’s disease: challenges and perspectives. Mol. Neurodegener.16, 76. 10.1186/s13024-021-00496-7
56
WangC.LiuH. (2022). Factors influencing degradation kinetics of mrnas and half-lives of micrornas, circrnas, lncrnas in blood in vitro using quantitative pcr. Sci. Rep.12, 7259. 10.1038/s41598-022-11339-w
57
WangL.ZhangL. (2020). Circulating exosomal miRNA as diagnostic biomarkers of neurodegenerative diseases. Front. Mol. Neurosci.13, 53. 10.3389/fnmol.2020.00053
58
WangZ. Y.SunM. H.ZhangQ.LiP. F.WangK.LiX. M. (2023). Advances in point-of-care testing of microRNAs based on portable instruments and visual detection. Biosensors13, 747. 10.3390/bios13070747
59
WoodsM. L.LeonM.Perez-CarrascoR.BarnesC. P. (2016). A statistical approach reveals designs for the Most robust stochastic gene oscillators. ACS Synth. Biol.5, 459–470. 10.1021/acssynbio.5b00179
60
ZhangC.TsoiR.WuF.YouL. (2016). Processing oscillatory signals by incoherent feedforward loops. PLOS Comput. Biol.12, e1005101. 10.1371/journal.pcbi.1005101
61
ZhangD. Y.SeeligG. (2011). Dynamic DNA nanotechnology using strand-displacement reactions. Nat. Chem.3, 103–113. 10.1038/nchem.957
62
ZhangD. Y.WinfreeE. (2009). Control of DNA strand displacement kinetics using toehold exchange. J. Am. Chem. Soc.131, 17303–17314. 10.1021/ja906987s
63
ZhangY.XuG.LianG.LuoF.XieQ.LinZ.et al (2020). Electrochemiluminescence biosensor for miRNA-21 based on toehold-mediated strand displacement amplification with Ru(phen)32+ loaded DNA nanoclews as signal tags. Biosens. Bioelectron.147, 111789. 10.1016/j.bios.2019.111789
64
ZiemssenT.DerfussT.de StefanoN.GiovannoniG.PalavraF.TomicD.et al (2016). Optimizing treatment success in multiple sclerosis. J. Neurology263, 1053–1065. 10.1007/s00415-015-7986-y
Summary
Keywords
feed-forward loops, toehold-mediated strand displacement, multi-objective optimisation, multiple sclerosis, miRNA, threshold detection, iGEM
Citation
Verkuijlen RJ and Smith RW (2025) A model-based design strategy to engineer miRNA-regulated detection systems. Front. Syst. Biol. 5:1601854. doi: 10.3389/fsysb.2025.1601854
Received
28 March 2025
Accepted
17 July 2025
Published
14 August 2025
Volume
5 - 2025
Edited by
George P. Patrinos, University of Patras, Greece
Reviewed by
Saptarshi Sinha, University of California, San Diego, CA, United States
Siyuan Wu, James Cook University, Australia
Updates

Check for updates
Copyright
© 2025 Verkuijlen and Smith.
This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Robert W. Smith, robert1.smith@wur.nl
Disclaimer
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.