{'fraction_non_alphanumeric': 0.021761812490684155, 'fraction_numerical': 0.012272072340636955, 'mean_word_length': 4.5802606043803715, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Amount of Escape Estimation Based on Bayesian and MCMC Approaches for RNA Interference

The amount of short interfering RNA (siRNA) escaping from the endosome has a significant impact on the efficiency of RNAi. In general, the initial injected amount of siRNAs during the experiment is known, and also the amount of siRNAs after the experiment can be revealed by the level of mRNA measured. However, it is impossible to measure the amount of siRNAs that escape from the endosome and really take part in the chemical reaction of RNAi by detecting the biological organism and its tissues. Inspired by the bottleneck effect in the virus, we introduce the Bayesian approach to infer the amount of escape based on a single type and multiple types of siRNA, respectively. With the consideration of the large calculation quantity of the accurate posterior distribution and the unavailable analytic expression of the likelihood function, our article proposes to take samples by the improved Markov chain Monte Carlo (MCMC) method. The article takes the silencing gene of the synthesis of chitin and the interfering multiple target oncogene as numerical examples to show that our improved MCMC method has higher operation efficiency compared to the Bayesian approach. Our research models siRNA endosome escape using statistical methods for the first time. It perhaps provides a theoretical basis to decrease the cost of a biotic experiment for the future and the standardized statistical approaches for the amount of escape estimation.


INTRODUCTION
RNAi refers to a highly conserved biological process that recognizes double-stranded RNA (dsRNA) in the cell to induce the specific degradation of homologous mRNA during evolution. 1 Endogenously expressed long dsRNA is first cleaved into short interfering RNA (siRNA) by the enzyme, such as Dicer, that is the component of a gene-silencing mechanism, and then the short RNA molecules are exploited as guides to target homologous RNA species. 2,3 The specific suppression of gene expression possibly actualizes through injecting or feeding with dsRNA. The introduction of siRNA into insect cells and silencing of target genes expression offer a new potential tool for the biological pest control method. 4 For example, the RNAi pathway could be applied to reduce the breeding of lepidopteran and coleopteran insect pests via restraining the planta expression, 5 and Mao et al. 6 provide a strategy to impair larval tolerance of gossypol by interfering a cotton bollworm RNA. As a highly efficient technology, RNAi has also developed rapidly in the field of infectious disease and tumor gene therapy, 7,8 and it can cure humans with various diseases that traditional drugs cannot, such as chronic hepatitis B virus. 9 In addition, individualized treatment schemes can be designed according to different conditions of patients.
The significant barrier for efficient siRNA uptake lies in the plasma membrane. In spite of the small size of siRNA molecules, they are still prevented from crossing biological membranes because of their negative charge and hydrophilicity. The procedure of the intracellular transportation of siRNAs begins with early endosomal vesicles. Subsequently, with the fusion of these early endosomes and sorting endosomes, siRNAs are transferred to the late endosomes. Only a small part of siRNAs could escape from the endosomes, and another part with the endosomal contents is removed to the lysosomes. The lysosomes that contain various nucleases acidify the endosomal content, and the siRNAs are degraded in turn. Figure 1 provides a schematic diagram that describes the process of the uptake and intracellular trafficking of a targeted siRNA. So, in order to avert lysosomal degradation, siRNAs have to escape from the endosomes and get into the cytosol, where they will associate with the RNAi mechanism. 10 Besides, it has been found that some of the generated siRNAs are not directly derived from the cleavage of dsRNA but rather, from a chain reaction of RNA polymerase. With the allowance of a single strand of siRNA as a primer and the target mRNA as a template, this reaction amplifies the target mRNA under the action of RNA-mediated RNA polymerase (RdRP) and generates a new siRNA subpopulation. 11 These, in turn, would continue to react to the target mRNA and degrade it. 12 This cyclical amplification process of RNAi explains the reason why a small amount of dsRNA can induce strong genesilencing effects.
We find that the process of siRNA delivery resembles the biological effect called bottleneck. The bottleneck describes the phenomenon that the number of individuals in a group is reduced drastically or even extinct due to drastic changes in the environment. When we inject a certain amount of siRNA into a pest, only a small fraction of the siRNA can across the plasma membrane and participate in the RNAi, and the remaining siRNAs will be degraded. The lower amount of escaping siRNA (commonly known as bottleneck size) will lead to a form of a new population by the amplification process. 13 Accurate quantification of the amount of escape for RNAi is vital for several reasons. First, the estimation of the amount of siRNAs escaping from the endosome helps us to research the biological mechanism of endosomal escape more definitively. Second, the knowledge of the amount of siRNAs of escape in RNAi processes is important to design rationally the strategies that optimize the amount of siRNA to interfere with the target RNA. Finally, the amount of escape impacts the levels of the types that can escape from the endosome into the cytosol when we inject multiple types of siRNA and thereby, impact the effect of interference.
Bottleneck has been extensively researched by many articles that mostly focus on the qualitative analysis of transmission bottleneck sizes, 14 and Abel et al. 15 provide a biologically motivated introduction to bottlenecks. Sobel et al. 16 use the deep-sequencing data to construct the likelihood expression of transmission bottleneck on the basis of the beta-binomial sampling method. Inspired by the above opinions with bottleneck, new ideas aiming at gauging the escaping amounts of siRNA for a single type and multiple types are suggested, respectively. After the observed data are simulated by the Gillespie algorithm, the probability distributions of escaping amounts of siRNA are estimated by means of two algorithms, consisting of the Bayesian approach and the nearest neighbor method. 17 However, both algorithms are inefficient in the course of actual implementation, because the multiple invoking and running of the Gillespie algorithm take much time. So we provide an alternative approach to sample the escaping amounts of siRNA based on the Markov chain Monte Carlo (MCMC) method and take the means of samples as the estimation of escaping amounts to improve the speed of the computer. Finally, comparisons indicate that the estimations inferred by both Bayesian and MCMC methods approximate the true value.

Silence Gene Controlling the Synthesis of Chitin
The oriental migratory locust is a crucial pest in agriculture. 18 Recently, the locust plague has broken out more frequently and severely in China. 19 As we know, the growth and development of locust strictly depend on the biosynthesis and degradation of chitin, which is absent in plants and vertebrates. So, chitin metabolism represents an attractive target for developing safe and effective insecticides. 20 RNAi can be used to silence genes that control the synthesis of chitin, sequentially leading to the death of locusts. After siRNAs are injected into the locust, they are governed by stochastic processes, including amplification, degradation, immigration, and emigration, which are dominated by a parameter set q = fa;l;m;sg. Let SðtÞ be the amount of the current siRNAs. Then, four stochastic processes are modeled by four biochemical reactions as follows: Next, the biological significance of the construction and parameters in Equation 1 are presented.
a is the amplification rate of siRNAs that have escaped. Equation 1a means that given the current amount SðtÞ, a unit of new siRNA is generated in the time interval ðt; t + dtÞ with probability aSðtÞdt. l is the degradation rate of siRNA due to the endocytosis. Equation 1b, represents that a unit of siRNA is degraded by lysosomes with probability lSðtÞdt in the time interval ðt; t + dtÞ for given the current states SðtÞ. m is the immigration rate of a new siRNA molecule. Equation 1c reveals that a unit of siRNA immigrates in our system from the neighboring cells with probability mdt in the time interval ðt;t + dtÞ. s is the emigration rate of siRNA. Equation 1d shows that siRNA will decrease one unit with the emigration of siRNA into the Take the parameter values a = 0:6; l = 0:3; m = 0:6; s = 0:23, for example, when the initial value is given by Sð0Þ = 5, simulations for the dynamic of the siRNA by the Gillespie algorithm are illustrated in Figure 2. So, the value at Dt = 12h could be recorded as our observation data s 2 being the amount of siRNA after amplification.
Next, the above observation data s 2 are employed to estimate the amount of escape s 1 or its posterior distribution pðs 1 j s 2 Þ and meanwhile, demonstrate the efficacy of Algorithm 1 and Algorithm 2 for the single type of siRNA.
For targets s Ã 1 = 7;s Ã 1 = 70, and s Ã 1 = 700, we obtained the posterior distributions pðs 1 j s 2 Þ of the escaping amount by Algorithm 1 in Figures 3A-3C. Furthermore, we take their modes 9, 67, and 687 as the estimations of the escaping amount, respectively. For the same targets, the samples of the escaping amount are displayed in Figures 4A-4C by Algorithm 2, and their means are estimated as 5, 78, and 687 after burn-in. Obviously, the two kinds of estimations fit the targets very well. This indicates that the two algorithms are efficient.

Interfere Multiple Target Oncogene
Related studies have found that the cancerization of normal cells is the consequence of interaction of multiple genes. However, conventional therapies, which are only targeted toward a single gene mostly, cannot completely inhibit the growth of tumors. It is obvious that RNAi technology can be utilized to silence gene. Yin et al. 22 suggested that injecting multiple types of siRNA can specifically interfere with multiple target oncogenes simultaneously and thereby inhibit the growth and proliferation of cancer cells synergistically.
Consequently, for multiple types, a hypothesis is given that we inject seven types of siRNA v ½1 0 ; v ½2 0 ; :::; v ½7 0 for gene therapy. Then, the observation data v 2 could be simulated by the Gillespie algorithm, as previously mentioned. Algorithm 3 and Algorithm 4 are applied to estimate the amount of escaping siRNAs and verify the efficacy of these two methods by the following steps.
Estimated posterior distributions pðs 1 j v 0 ; v 2 Þ by Algorithm 3 are shown in Figures 5A-5C for targets s Ã 1 = 7; s Ã 1 = 70, and s Ã 1 = 140. The modes, as the estimations of the escaping amount, are 9, 65, and 135, respectively. For the same targets, we perform 10,000 samples by Algorithm 4 and report in Figures 6A-6C. After burn-in, we get the estimations 9, 69, and 133 by calculating their means. It can be seen that our predicted results approximate accurately to real ones.

DISCUSSION
The amount of siRNAs escaping from the endosome is one of the important essentials dominating the efficiency of RNAi, but it is intractable to be observed and calculated in experiments. In this paper, two methods are proposed to estimate the amount of escape in terms of the knowledge of the dynamics during amplification from the amount after the reaction and the amount of injection. One is to estimate the posterior distribution of escaping the amount according to the Bayesian approach; the other one is to get the samples of the escaping amount by the MCMC method and to use the mean of samples as an estimate. For the traditional Bayesian  approach, we present the specific algorithms combined with the nearest neighbor method, which is used for the estimation of pðs 2 j s 1 Þ. For the MCMC method, the acceptance probability of the Metropolis-Hastings (MH) algorithm is controlled by the distance function between the simulation with the observed data. Furthermore, with the contraposition of the single type of siRNAs and multiple types of siRNAs, the algorithms of the estimate of the escaping amount are given, respectively. To inspect the validity of our algorithms, two examples on the silencing gene for the synthesis of chitin and blocking multiple target oncogenes are derived. Our pursuit offers statistical ways to infer the exact amount of siRNAs participating in the actual RNAi reaction. Meanwhile, it perhaps provides a theoretical basis to decrease the cost of the biotic experiment for the future.
Even so, there are still some problems worth exploring further. First, the MCMC method failed to estimate the posterior distribution that could express the uncertainty through the variance of the distributions, although it improves the efficiency. It indicates that a more comprehensive method that takes into account the accuracy of estimation, efficiency, and expression of uncertainty together is required. Besides, the estimation of the bottleneck size is only built on the assumption that the dynamics during amplification are known. When the partial data are missing, how to estimate the amount of escape and the parameters together is the problem for further consideration. In future research, we will try to find the solutions to these problems.

Single Type of siRNA
In general, we only introduce a single type of siRNA aimed at a specific RNA into the organisms. The processes for which siRNAs escape from the endosome and amplify intracellularly have been described in the first part, and now, we picture them in Figure 7. Define the initial injected amount of siRNAs as s 0 , the amount of siRNAs that escape from endosome as s 1 , and the amount of siRNAs after amplification as s 2 (Figure 7). Obviously, s 1 %s 2 . Then, on the premise of the amount of siRNAs after amplification, Bayesian inference or MCMC can be applied to estimate the posterior distribution of the escaping amount of siRNAs, as well as their value.

Bayesian Inference
According to the Bayesian framework, the amount of siRNAs escaping from the endosome can be estimated by the posterior probability distributions. Given the observations of the amount after amplification, the distribution is given by pðamount of escape ðs 1 Þ j amount after amplification ðs 2 ÞÞ: The merit of the use of the Bayesian approach is that we not only could get the estimates of the most probable amount of escape (in terms of the modes of the distribution), but also, we could be aware of the uncertainty via the variance of the distributions. Then, the posterior probability pðs 1 j s 2 Þ is given by  Then, the posterior distribution pðs 1 j s 2 Þ can be obtained through estimating all of the probability pðs 1 j s 2 Þ for s 1 = 1; :::; s max , where s max is the maximum of escaping amount s 1 . The detailed process is shown as follows.
First, starting from s 1 , we perform n simulations using the Gillespie stochastic algorithm, 21 according to a parameter set q for the dynamics, and obtain the finite simulating samples of s 2 after time interval Dt: È ðs 2 Þ 1 ; ðs 2 Þ 2 ; :::; ðs 2 Þ n É ; from which pðs 2 j s 1 ; qÞ is estimated using the nearest neighbor method, 17 which is a classical nonparametric estimation method.
Second, with the substitution of all probabilities pðs 2 j ,; ,Þ into Equation 3, one gets the estimation of the probability distribution pðs 1 j s 2 Þ .
In detail, the algorithm for estimating distribution pðs 1 j s 2 Þ is given as follows.
Algorithm 1 implies that the Gillespie algorithm runs n times when the loop executes one time. It reveals that Algorithm 1 is time consuming if simulating time n is large. So, in order to improve the running efficiency of program, we adopt the MCMC method to estimate the escaping amount of siRNAs.

MCMC Method
The MCMC method includes Gibbs and MH, which are techniques simulating the random variables by using the Markov chain. 23 In this paper, we choose MH to sample single variable s 1 , rather than Gibbs from the target distribution, being the conditional distribution of interest. Here, the target distribution, that is, posterior distribution pðs 1 j s 2 Þ in Equation 2, is proportional to the product of prior pðs 1 Þ and likelihood pðs 2 j s 1 Þ.
From the ideas of MCMC, we need to compute the acceptance probability aðs ' Output: the probability pðs 1 j v 0 ; v 2 Þ when s 1 = 1; :::; s max .
Sample v 1 from the multivariate hypergeometric distribution with s 1 , v 0 4.
Simulate  Now, the procedure of sampling s 1 by MCMC methods is listed as follows.

Multiple Types of siRNA
With the consideration of injecting multiple types of siRNAs to affect different target RNAs, the stochastic process of siRNAs is shown in Figure 8. Assume that we inject m types of siRNA for which the initial injected amount consists of v Therefore, the incorporation of Equations 8 and 9 yields Assume that all types of siRNAs are phenotypically identical and have the same probability of escaping from the endosome. Then, the distribution v of s 1 could be considered as sampling randomly without replacement from the initial injected amount with distribution v 0 . So, we can select the amount of escape from a multivariate hypergeometric distribution with v 0 and sumðv 1 Þ = s 1 . The probability of drawing v 1 from v 0 is given by In reality, components of v 2 are simulated by the Gillespie algorithm in view of parameter vector q. So, for convenience, pðv 2 j v 1 Þ is denoted by pðv 2 j v 1 ; qÞ, which is factorized in accordance with the independence between each type of siRNA as follows: Then, pðv ½i 2 v ½i 1 ; qÞ could be estimated the same way that we estimate pðs 2 j s 1 ; qÞ, used in Algorithm 1.
The acquisition of pðv 1 j v 0 Þ and pðv 2 j v 1 Þ that are desired for Equation 10 has been solved in the previous segment, but we should count all of the summands when v 2 gets every possible value, such that One key problem is that all possible values of v 1 grow superexponentially with s 1 when we give a value of s 1 . 26 Now, we face a combinatorial and computational challenge, and so a replaceable approach is required.
To avoid the combinatorial problem, the more probable configuration of v 1 , such as the modes of v 1 , could replace the summands that consider all possibilities of v 1 in Equation 10. Requena et al. 27 have elaborated an algorithm to solve this question, but now, we provide a simpler sampling method that is to sample points v 1 randomly from multivariate hypergeometric distribution pðv 1 j v 0 Þ for enough times so that most of these points would be adjacent to the modes. The concrete execution of the sampling procedure is shown in Algorithm 3.
Likewise, as discussed in the context above, there are problems of efficiency with this approach. Therefore, it is tempting to attempt to use the MCMC method.

MCMC Method
Multiple types are also appropriate for the MCMC method. Similar to the single type, our target distribution is posterior distribution pðs 1 j v 0 ; v 2 Þ now. From Equation 9, we get In order to go to the acceptance probability, first, we should draw v ' 1 from the multivariate hypergeometric distribution with s 1 and given v 0 . Afterward, simulate one vector of v ' 2 from v ' 1 after Dt by the Gillespie algorithm, and then, the distance between the given value v 2  The exact process of the MCMC method is described in Algorithm 4.

AUTHOR CONTRIBUTIONS
Y.P. conceived the project and designed the frame of this paper; T.L. and C.L. finished mathematical analyses, performed simulations and wrote the first draft; M.Y. polished, revised the last draft. All authors contributed to the manuscript and approved the final manuscript.

CONFLICTS OF INTEREST
The authors declare no competing interests.