Noise propagation with interlinked feed-forward pathways

Functionally similar pathways are often seen in biological systems, forming feed-forward controls. The robustness in network motifs such as feed-forward loops (FFLs) has been reported previously. In this work, we studied noise propagation in a development network that has multiple interlinked FFLs. A FFL has the potential of asymmetric noise-filtering (i.e., it works at either the “ON” or the “OFF” state in the target gene). With multiple, interlinked FFLs, we show that the propagated noises are largely filtered regardless of the states in the input genes. The noise-filtering property of an interlinked FFL can be largely derived from that of the individual FFLs, and with interlinked FFLs, it is possible to filter noises in both “ON” and “OFF” states in the output. We demonstrated the noise filtering effect in the developmental regulatory network of Caenorhabditis elegans that controls the timing of distal tip cell (DTC) migration. The roles of positive feedback loops involving blmp-1 and the degradation regulation of DRE-1 also studied. Our analyses allow for better inference from network structures to noise-filtering properties, and provide insights into the mechanisms behind the precise DTC migration controls in space and time.

no activation in UNC-5, which corresponds to the retarded phenotype. Mutant blmp-1;daf-12 has an oscillatory activity for UNC-5, which we further study with stochastic simulation.
The phenotype distribution from stochastic simulation, for the genotypes that are not included in the maintext, is included in Fig. S5. We observe that DTC turning distribution in lin-29 and daf-12 single mutants is slightly wider compared to WT (Fig. 7c in main text). This result correlates with the higher noise we observed in these mutants at the time of DTC turning (late L3 stage) (Fig 6a in main text). Single mutation in blmp-1 results in a precocious distribution where DTC turning observed at early L3 stage and with some at late L2. The three double mutants, dre-1;daf-12, lin-29;daf-12, and lin-29;dre-1 do not show UNC-5 expression even at 20 hr.
In Table S1 we list the results of testing over the logic combination for the regulation of blmp-1 from the three upstream genes and an auto-activation. Reported are the fractions of parameter sets that passed the screening of wild type or those also passed the 6 mutants as described in the Method section of the main text. It is seen that logic No.2 has the highest rates in both cases, and therefore, it is used for all the computational analysis in the present work.

Alternative random screening for parameters
We have tested with a random search (Monte Carlo simulation) with independent production and degradation rates for each gene (Table S2), and it was repeated with two different burst sizes. Again we screen for parameter sets that can reproduce the wild type using deterministic ODE by requiring a clear distinction between the ON and OFF state level as described in Method section(page 15) of the main text.
We confirmed that the noise filtering capacity of subnetworks through IFFL and PFL, built from this more general random sampling of parameters, is similar to the conclusions derived from the complete scanning performed and reported in the main text. From Fig.S6, S7 and S8, it is seen that the changes from case to case remain the same, but the observed FF levels are higher with a wider distribution, simply due to the different steady state obtained from various production rates sampled. Similar results are also observed (Fig.S9, S10 and S11) when the burst size was set at 64.

Additional remarks on parameter settings
The time scales of changes are mainly determined by the degradation time (γ's). 1, 2 In a dynamic system for developments such as the currently studied DTC cell migration, gene expression changes at various levels through the development stages, measured in hours. Therefore, the degradation time should allow such change to take place with a time scale of hours, which set a lower bound for γ's. Very large γ's allow fast switch for the genes, which would allow the same conclusion in modeling the simulation, but very fast switches would make numerical propagation of the ODE's more time-consuming but with essentially the same results, and therefore we pick a range of γ's according to these constraints. BLMP-1 has additional degradation effect by DRE-1 as observed in experiments 3, 4 supporting that BLMP-1 self-degradation rate can be lower than other components in the network. So in all simulations in the present work, we assume BLMP-1 basal degradation is slower than other regulating components such that DRE-1 degradation is effective in simulation and the experimental observations can be reproduced. Again our results showed that the main conclusion on noise-filtering is unaffected as long as it is set to be slower than other components.
The level of gene expression (k p /γ) has been quantified for many different systems. 5, 6 While it is not known for the specific development state that we focused on, we can simply allow it to very across a reasonable range as set by experimental results 5,6 in this Monte Carlo sampling. However, we note that none-dimensionalization allows us to choose arbitrary units for each component for deterministic simulation, as people have published excellent works with this non-dimensional setting. 7,8 Therefore, our deterministic results in the main text is not affected by the setting in the gene expression. The new results from Monte Carlo sampling indicate that the noise-filtering behavior is indeed not affected.
The burst size b has a direct consequence on the size of fluctuation in each component. And the fluctuation is relative to gene expression level, and thus, the gene expression level does have an effect on the noise amplitude. In the main text, we simplified the simulation by using a fixed overall expression level and a fixed burst size. In the supplementary we included new Monte Carlo sampled parameter search with a large variation in the gene expression rates but with a fixed burst size. We choose the latter setting to test simply because the effects of varying both is mainly in the variation of fluctuation (Fano factor) of each gene, which can be equally achieved by varying the expression rates only. The fact that both simulation leads to similar conclusions in noise-filtering indicates that the level of noise (and factors behind such noises) is largely independent to the noise-filtering mechanism. We note that, in an additional test, a change in the burst size (from 16 to 64) leads to the same overall dynamics, the same changes in the fluctuation, and the same conclusions in noise filtering, which further confirms our conclusion here.
The threshold value (K's) of each gene regulation has the same dimension as the corresponding regulator's abundance. For a gene regulation to be effective, it has to be turned ON and OFF as observed and concluded in experiments. Therefore, the  threshold values has been chosen according to the ON and OFF-state level of the regulators in both tests. We note that a wider choice of the threshold values is not that necessary because it would simply create an all-ON (for very small K's) or all-OFF (for very large K's) situation through out the dynamics of development, and such situation means a change in regulation is not possible, which is often not what experimentally observed, and can be ruled out safely.
In the parameter scanning presented in the main text, the combination of all 13 parameters yields 1.48 × 10 10 different combination already. In that simulation, we obtained 3.08 × 10 7 combinations that pass the wild-type requirements, and a randomly chosen 1000 sets were used for the current work, which offers a very general ground for the conclusions. In the new round of Monte Carlo simulation (with b equals 16), the number of independent parameters was increased to 23, with a broader parameter setting, the chance of finding a suitable parameter set dropped from 0.2 % to 0.006 %. Again the same conclusion can be drawn for this new setting.