Messenger Ribonucleic Acid Sequence Design Method, Device, Computing Device and Storage Medium
Optimizing mRNA sequences through growth method and cluster cluster search method has solved the problem of failure to take into account the stability and expression of the full-length and CDS regions in the prior art, and improved the stability and expressability of mRNA sequences, reducing storage and transportation costs.
Patent Information
- Application Number
- CN202311113680.3
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- Filing Date
- 2023-08-31
- Publication Date
- 2025-07-29
- Estimated Expiration
- 2043-08-31
AI Technical Summary
The existing mRNA sequence design methods fail to effectively take into account the minimum free energy of the full length and the codon usage efficiency of the CDS region, resulting in the instability of synthetic mRNA molecules, which requires low temperature storage and increased storage and transportation costs.
Codons were added one by one by growth method, combined with the calculation and scoring and sorting of codon usage efficiency of the CDS region codon, the elimination method was used to exclude bad sequences, and representative sequences were selected through the cluster cluster search method to optimize mRNA sequences in 5’UTR and 3’UTR environments.
The stability and expressibility of mRNA sequences are achieved while optimizing, reducing storage and transportation requirements, improving the diversity and innovation of design sequences, and reducing the computational volume.
Smart Images

Figure CN119541638B_ABST
Abstract
Description
Technical Field
[0001] The present invention relates to the technical field of biopharmaceuticals, and particularly relates to a messenger ribonucleic acid sequence design method, device, computing device, and storage medium. Background Art
[0002] Messenger ribonucleic acid (mRNA) is a single-stranded molecule responsible for transmitting genetic information from deoxyribonucleic acid (DNA) to ribosomes, which then decode the genetic information and synthesize proteins. In recent years, the technical means of synthesizing mRNA via in vitro transcription (IVT) technology has gradually matured. Biomedical companies use the mRNA technology platform and can deliver mRNA into the human body using delivery vectors with targeting properties, high efficiency, and low immunogenicity. Human cells utilize the synthesized mRNA to express specific proteins. Since the synthesized mRNA is structurally similar to natural mRNA, patients can produce therapeutic proteins in their own bodies, thus reducing various troubles caused during the complex manufacturing process of recombinant proteins. Due to the diversity of mRNA-based therapies, new drug developers use synthetic mRNA as a disease treatment tool and have made progress in the research and development of cancer immunotherapy, gene therapy, and the control of infectious diseases. The mRNA technology has obvious advantages in terms of efficiency and protective power in vaccine research and development. Both Moderna and BioNTech developed and launched mRNA vaccines within a year, creating history.
[0003] When designing mRNA sequences, biomedical companies often do not directly adopt wild-type mRNA sequences but optimize the design of mRNA sequences according to the requirements of activity and drug-likeness. Traditional mRNA sequence design methods usually only consider codon optimization, that is, using codons with a higher frequency of use in human cells to optimize the CDS sequence encoding proteins, so that the synthesized mRNA sequence has higher expressibility compared to the wild-type mRNA sequence. However, this design method that only considers the codon usage frequency lacks consideration of physicochemical properties, making the designed synthetic mRNA molecules unstable and prone to degradation in vitro. Therefore, low-temperature (-20 / -80°C) storage conditions are required, increasing the storage and transportation costs.
[0004] To address the above problems, Huang Liang et al. from Baidu Research Institute in the United States proposed a new algorithm for mRNA sequence design: Linear Design in 2020. In this algorithm, the optimization of the CDS region of synthetic mRNA is set as a multi-objective optimization problem: that is, simultaneously optimizing the minimum free energy (MFE) and codon adaptation index (CAI) of the CDS region. By defining a scoring function S = MFE CDS + λ|p|(-logCAI CDS ), Linear Design simultaneously optimizes the minimum free energy and codon adaptation efficiency of the CDS region of an mRNA sequence to obtain a synthetic mRNA molecule with both good expressivity and stability. However, in addition to the CDS region encoding proteins, synthetic mRNA molecules for pharmaceutical industrial applications also include non-coding regions such as 5'UTR and 3'UTR, and these non-coding regions also have important effects on the expression and stability of mRNA molecules.
[0005] Therefore, there is an urgent need for a method for designing messenger ribonucleic acid sequences that can simultaneously optimize the minimum free energy of the full length of an mRNA sequence and the codon adaptation efficiency of the CDS region. Summary of the Invention
[0006] The purpose of the present invention is to provide a method, device, computing device, and storage medium for designing messenger ribonucleic acid sequences, which can simultaneously optimize the minimum free energy of the full length of an mRNA sequence and the codon adaptation efficiency of the CDS region, ensuring the stability and expressibility of the mRNA sequence.
[0007] To solve the above technical problems, an embodiment of the present invention discloses a method for designing messenger ribonucleic acid sequences, including the following steps:
[0008] Step S1: Add 1 codon from all possible codons corresponding to the first amino acid of the encoded protein to the end of the given messenger ribonucleic acid sequence respectively to generate a group of initial sequences, and put the initial sequences into a sequence pool; where only 1 codon is added to the end of the sequence, and each codon is only added once;
[0009] Step S2: Determine whether the number of sequences in the sequence pool is greater than the search depth D;
[0010] If not, execute step S3: Add one codon from all possible codons corresponding to the next amino acid of the encoded protein to the end of each sequence in the sequence pool, thereby generating a new set of sequences, and return to execute step S2; wherein, only one codon is added to the end of the sequence, and each codon is only added once;
[0011] If so, execute step S4: Calculate the full-length minimum free energy MFE of each sequence in the sequence pool whole and the codon adaptation index CAI of the CDS region CDS ;
[0012] Step S5: Calculate the score of each sequence in the sequence pool according to the full-length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS , calculate the scores of each sequence in the sequence pool, and sort all the sequences in the sequence pool according to the scores;
[0013] Step S6: Filter out the sequences that do not meet the requirements in the sequence pool or eliminate the sequences ranked at the end in the sequence pool according to the preset filtering rules or elimination rate, and the remaining sequences are retained;
[0014] Step S7: Cluster the retained sequences into D categories according to the search depth D, and retain the sequence with the highest score in each category in the sequence pool, and eliminate the remaining sequences;
[0015] Step S8: Determine whether all the codons corresponding to all the amino acids of the encoded protein have been added;
[0016] If not, return to execute step S3; if so, execute step S9: Add the 3'UTR sequence to the end of each sequence in the sequence pool;
[0017] Step S10: Calculate the full-length minimum free energy MFE of each sequence in the sequence pool after adding the 3'UTR sequence whole and the codon adaptation index CAI of the CDS region CDS ;
[0018] Step S11: Calculate the score of each sequence in the sequence pool after adding the 3'UTR sequence according to the full-length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS , calculate the scores of each sequence in the sequence pool after adding the 3'UTR sequence, and sort all the sequences in the sequence pool after adding the 3'UTR sequence according to the scores;
[0019] Step S12: According to the number R of the selected sequences, cluster all the sequences in the sequence pool after adding the 3'UTR sequences into R categories, and recommend the sequence with the highest score in each category as the designed messenger ribonucleic acid sequence.
[0020] An embodiment of the present invention also discloses a ribonucleic acid sequence design device, including:
[0021] A first sequence generation unit, configured to add one codon from all possible codons corresponding to the first amino acid of the encoded protein to the end of the given messenger ribonucleic acid sequence respectively to generate a set of initial sequences, and put the initial sequences into a sequence pool; wherein, only one codon is added to the end of the sequence, and each codon is added only once;
[0022] A first judgment unit, configured to judge whether the number of sequences in the sequence pool is greater than the search depth D;
[0023] A second sequence generation unit, configured to, when the first judgment unit judges that the number of sequences in the sequence pool is not greater than the search depth D, add one codon from all possible codons corresponding to the next amino acid of the encoded protein to the end of all the sequences in the sequence pool respectively, so as to generate a set of new sequences; wherein, only one codon is added to the end of the sequence, and each codon is added only once;
[0024] A first calculation unit, configured to, when the first judgment unit judges that the number of sequences in the sequence pool is greater than the search depth D, calculate the full-length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS ;
[0025] A second calculation unit, configured to calculate the score of each sequence in the sequence pool according to the full-length minimum free energy MFE whole output by the first calculation unit CDS and the codon adaptation index CAI of the CDS region, and sort all the sequences in the sequence pool according to the scores;
[0026] A filtering unit, configured to filter the sequences that do not meet the requirements in the sequence pool or eliminate the sequences ranked at the back in the sequence pool according to a preset filtering rule or elimination rate, and the remaining sequences are retained;
[0027] A first clustering unit, configured to cluster the retained sequences into D categories according to the search depth D, and retain the sequence with the highest score in each category in the sequence pool, and eliminate the remaining sequences;
[0028] A second judgment unit that judges whether all the codons corresponding to the amino acids of the encoded protein have been added;
[0029] The second sequence generation unit is further configured to, when the second judgment unit judges that not all the codons corresponding to the amino acids of the encoded protein have been added, add one codon from all the possible codons corresponding to the next amino acid of the encoded protein to the end of each sequence in the sequence pool, so as to generate a set of new sequences; wherein, only one codon is added to the end of the sequence, and each codon is only added once;
[0030] A 3' UTR sequence addition unit that, when the second judgment unit judges that all the codons corresponding to the amino acids of the encoded protein have been added, adds a 3' UTR sequence to the end of each sequence in the sequence pool;
[0031] A third calculation unit that calculates the full-length minimum free energy MFE of each sequence in the sequence pool after adding the 3' UTR sequence whole and the codon adaptation index CAI of the CDS region CDS ;
[0032] A fourth calculation unit that, based on the full-length minimum free energy MFE output by the third calculation unit whole and the codon adaptation index CAI of the CDS region CDS , calculates the score of each sequence in the sequence pool after adding the 3' UTR sequence, and sorts all the sequences in the sequence pool after adding the 3' UTR sequence according to the score;
[0033] A second clustering unit that selects the number R of sequences, clusters all the sequences in the sequence pool after adding the 3' UTR sequence into R categories, and recommends the sequence with the highest score in each category as the designed messenger ribonucleic acid sequence.
[0034] An embodiment of the present invention also discloses a messenger ribonucleic acid sequence design calculation device, including:
[0035] A memory for storing computer-executable instructions; and,
[0036] A processor that, when executing the computer-executable instructions, implements the steps in the above method.
[0037] An embodiment of the invention also discloses a computer-readable storage medium, in which computer-executable instructions are stored, and when the computer-executable instructions are executed by a processor, the steps in the above method are implemented.
[0038] Compared with the prior art, the main differences and effects of the embodiments of the present invention are as follows:
[0039] It can simultaneously optimize the minimum free energy of the full length of an mRNA sequence and the codon usage efficiency of the CDS region, ensuring the stability and expressibility of the mRNA sequence.
[0040] Furthermore, by using the growth method, the computational amount is effectively controlled, avoiding the exhaustion of computing power.
[0041] Furthermore, by using the elimination method, defective sequences or sequences with poor scores can be avoided, ensuring the quality of the designed sequences.
[0042] Furthermore, through cluster search and recommendation, the diversity of the designed sequences is improved, and the innovation of the designed sequences is enhanced.
[0043] A large number of technical features are recorded in the specification of this application, distributed in various technical solutions. If all possible combinations of technical features (i.e., technical solutions) of this application are to be listed, the specification will be too long. To avoid this problem, each technical feature disclosed in the above-mentioned invention content of this application, each technical feature disclosed in the following embodiments and examples, and each technical feature disclosed in the drawings can be freely combined with each other to form various new technical solutions (these technical solutions are all regarded as having been recorded in this specification), unless the combination of such technical features is technically infeasible. For example, in one example, features A + B + C are disclosed, and in another example, features A + B + D + E are disclosed, and features C and D are equivalent technical means that play the same role. Only one of them can be used technically and it is impossible to use both at the same time. Feature E can be combined with feature C technically. Then, the solution of A + B + C + D should not be regarded as having been recorded due to technical infeasibility, while the solution of A + B + C + E should be regarded as having been recorded. Brief Description of the Drawings
[0044] Figure 1 is a schematic flowchart of a method for designing a messenger ribonucleic acid sequence in the first embodiment of the present invention;
[0045] Figure 2 is a comparison schematic diagram between the first embodiment of the present invention and the prior art;
[0046] Figure 3 is a schematic flowchart of a method for designing a messenger ribonucleic acid sequence in a preferred embodiment of the present invention;
[0047] Figure 4 is for optimizing the objective function MFE whole + K * CAI CDSSchematic diagram of mRNA design of Spike protein under specific UTR constraints using the Beam Search method;
[0048] Figure 5 is under the optimized objective function MFE whole +K*CAI CDS Schematic diagram of mRNA design of Spike protein under specific UTR constraints using the Beam Cluster Search method;
[0049] Figure 6 Schematic diagram of the comparison of the identity of the CDS regions of the Spike protein mRNA sequences designed using the Beam Search and Beam Cluster Search methods;
[0050] Figure 7 Schematic diagram of the comparison of the identity of the CDS regions of the Spike protein mRNA sequence (K = 4, D = 100, drop rate = 0.02) designed by the Beam Cluster Search method with those of the marketed vaccines of BioNTech and Moderna;
[0051] Figure 8 Schematic diagram of the comparison of the identity of the CDS regions of the Spike protein mRNA sequence (K = 4, D = 100) designed by the Beam Search method with those of the marketed vaccines of BioNTech and Moderna;
[0052] Figure 9 Schematic diagram of the influence of the Drop Rate (DR) on the diversity of the designed sequences;
[0053] Figure 10 Schematic diagram of the influence of the Depth (D) on the diversity of the designed sequences;
[0054] Figure 11 Schematic diagram of the structure of a messenger ribonucleic acid sequence design device in the second embodiment of the present invention. Detailed implementation mode
[0055] In the following description, many technical details are presented to enable the reader to better understand the present application. However, those of ordinary skill in the art can understand that the technical solutions claimed in the various claims of the present application can be implemented even without these technical details and various changes and modifications based on the following embodiments.
[0056] Explanation of some concepts:
[0057] mRNA: Messenger Ribonucleic Acid
[0058] MFE (minimum free energy): Minimum Free Energy
[0059] CAI (codon adaption index): Codon Adaptation Index
[0060] diversity: Diversity
[0061] To make the objectives, technical solutions and advantages of the present invention clearer, the embodiments of the present invention will be further described in detail below with reference to the accompanying drawings.
[0062] The first embodiment of the present invention relates to a method for designing a messenger ribonucleic acid sequence. Figure 1 It is a schematic flowchart of the method for designing the messenger ribonucleic acid sequence.
[0063] For a given protein sequence, the number of theoretically corresponding codon permutation and combination spaces is extremely large. To avoid exhausting computing power by enumerating all possibilities, we use the Growing Algorithm to search for optimized sequences.
[0064] Specifically, as Figure 1 shown, the method for designing the messenger ribonucleic acid sequence includes the following steps:
[0065] In step S1, add 1 codon from all possible codons corresponding to the first amino acid of the encoded protein to the end of the given messenger ribonucleic acid sequence respectively to generate a group of initial sequences, and put the initial sequences into a sequence pool; wherein, only 1 codon is added to the end of the sequence, and each codon is only added once.
[0066] In this embodiment, preferably, the given messenger ribonucleic acid sequence is a 5' UTR sequence.
[0067] That is to say, in this step S1, the starting sequence is a 5' UTR sequence.
[0068] Taking this 5' UTR sequence as the starting input, and then based on this 5' UTR sequence, add all possible codons corresponding to the first amino acid of the encoded protein to the end of the 5' UTR sequence, thereby generating a group of initial sequences, and put the initial sequences into the sequence pool.
[0069] It should be noted that only 1 codon is added to the end of the 5' UTR sequence each time, and each codon among all possible codons corresponding to the first amino acid of the encoded protein is only added once.
[0070] In this embodiment, preferably, the number of all possible codons corresponding to one amino acid of the encoded protein is 1 - 6.
[0071] Correspondingly, in step S1, each time 1 codon is added to the end of the 5'UTR start sequence respectively. After each codon among all possible codons corresponding to the first amino acid of the encoded protein is added respectively, this 1 5'UTR start sequence grows into 1 - 6 initial sequences.
[0072] Thereafter, step S2 is entered: Determine whether the number of sequences in the sequence pool is greater than the search depth D.
[0073] If not, then step S3 is executed: Add 1 codon among all possible codons corresponding to the next amino acid of the encoded protein to the end of all sequences in the sequence pool respectively, so as to generate a group of new sequences, and return to execute step S2; wherein, only 1 codon is added to the end of the sequence, and each codon is only added once.
[0074] If so, then step S4 is executed: Calculate the full - length minimum free energy MFE of each sequence in the sequence pool whole and the codon adaptation index CAI of the CDS region CDS .
[0075] Thereafter, step S5 is entered: According to the full - length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS , calculate the score of each sequence in the sequence pool, and sort all sequences in the sequence pool according to the score.
[0076] In this embodiment, preferably, the score of the sequence is the weighted sum of the full - length minimum free energy MFE of the whole sequence whole and the codon adaptation index CAI of the CDS region CDs . The formula for the score of the sequence is expressed as:
[0077] Score = MFE whole + K * CAI CDS
[0078] wherein, Score is the score of the sequence, and K is the weight coefficient.
[0079] Furthermore, preferably, in step S5, the formula for the score of the sequence is further expressed as:
[0080]
[0081] Wherein, Score is the score of the sequence, K is the weight coefficient, m is the number of codons in the CDS region, and i represents the i-th codon in the CDS region.
[0082] In this step S5, according to the preset target optimized weight coefficient K value, calculate the score of each sequence, and sort them from low to high.
[0083] The value range of K is [0, +∞].
[0084] K = 0 means that the score of the above sequence only considers the minimum free energy of the full length of the sequence; K = +∞ means that the score of the above sequence only considers the codon usage efficiency of the CDS region; if K is a positive number, for example, K = 10.0, it means that the score of the above sequence should reach a balance between the minimum free energy of the full length of the overall sequence and the codon usage efficiency of the CDS region.
[0085] Then enter step S6: According to the preset filtering rules or dropout rate, filter out the sequences that do not meet the requirements in the sequence pool or eliminate the sequences ranked at the end in the sequence pool, and the remaining sequences are retained.
[0086] In this step S6, use the elimination method to exclude the sequences containing restriction sites / undesirable properties / or poor scores.
[0087] Then enter step S7: According to the search depth D, cluster the retained sequences into D categories, and retain the sequence with the highest score in each category in the sequence pool, and eliminate the remaining sequences.
[0088] In this embodiment, preferably, the search depth D can be any one of 50, 100, 200, and 500.
[0089] In this embodiment, take the search depth D as 100 as an example for illustration.
[0090] Of course, this is only a preferred embodiment of the present application. D can be any one of 50, 100, 200, and 500, or even higher values, and is not limited thereto.
[0091] The higher the D value, the more sequences enter the next round, and the better the diversity of the designed sequences. Of course, the computational cost is also greater.
[0092] In this step S7, among the retained sequences, cluster them according to the search depth. If D = 100, cluster them into 100 clusters, and then select the sequence with the highest score from each cluster and retain it in the sequence pool, and eliminate the remaining sequences.
[0093] Then, step S8 is entered: Determine whether all the codons corresponding to the amino acids of the encoded protein have been added.
[0094] If not, return to step S3 for execution.
[0095] And so on, repeat steps S3 - S8 until all the codons of the encoded protein have been added.
[0096] After multiple rounds of growth, the number of sequences retained in the sequence pool from the previous round is fixed at D. After each codon is grown, the number of sequences in the sequence pool becomes D - 6D. Therefore, D sequences must be selected to enter the next round of cycle.
[0097] If so, execute step S9: Add 3'UTR sequences to the ends of all the sequences in the sequence pool.
[0098] Finally, add the 3'UTR sequence, thus growing into a full - length mRNA sequence (5'UTR + CDS + 3'UTR). This full - length mRNA sequence includes not only the CDS region encoding the protein but also non - coding regions such as 5'UTR and 3'UTR, and these non - coding regions also have important effects on the expression and stability of the mRNA molecule.
[0099] Then, step S10 is entered: Calculate the full - length minimum free energy MFE of each sequence in the sequence pool after adding the 3'UTR sequence whole and the codon adaptation index CAI of the CDS region CDS .
[0100] Then, step S11 is entered: According to the full - length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS , calculate the score of each sequence in the sequence pool after adding the 3'UTR sequence, and sort all the sequences in the sequence pool after adding the 3'UTR sequence according to the score.
[0101] Similarly, preferably, in this step S11, the score of the sequence is the weighted sum of the full - length minimum free energy MFE whole of the whole sequence and the codon adaptation index CAI CDS of the CDS region. The formula for the score of the sequence is expressed as:
[0102] Score = MFE whole + K * CAI CDS
[0103] where Score is the score of the sequence and K is the weight coefficient.
[0104] Further, preferably, the formula for the score of the sequence is further expressed as:
[0105]
[0106] Wherein, Score is the score of the sequence, K is the weight coefficient, m is the number of codons in the CDS region, and i represents the i-th codon in the CDS region.
[0107] In this step S11, according to the preset target optimization weight coefficient K value, calculate the score of each sequence and sort them from low to high.
[0108] The value range of K is [0, +∞].
[0109] K = 0 means that the score of the above sequence only considers the minimum free energy of the full length of the sequence; K = +∞ means that the score of the above sequence only considers the codon usage efficiency of the CDS region; if K is a positive number, for example, K = 10.0, it means that the score of the above sequence should achieve a balance between the minimum free energy of the full length of the overall sequence and the codon usage efficiency of the CDS region.
[0110] Thereafter, enter step S12: According to the number R of the selected sequences, cluster all the sequences after adding the 3'UTR sequence in the sequence pool into R categories, and recommend the sequence with the highest score in each category as the designed messenger ribonucleic acid sequence.
[0111] Thereafter, end this process.
[0112] In summary, the embodiment of the present application proposes a method for designing a messenger ribonucleic acid (mRNA) sequence in an environment with 5'UTR and 3'UTR, which can simultaneously optimize the minimum free energy of the full length of an mRNA sequence and the codon usage efficiency of the CDS region, ensuring the stability and expressibility of the mRNA sequence.
[0113] In order to better understand the technical solution of this specification, a preferred embodiment is described below. The details listed in this preferred embodiment are mainly for easy understanding and do not limit the protection scope of the present application.
[0114] The technical solution of this preferred embodiment proposes a new method for optimizing the design of mRNA sequences. Different from the prior art, which only simultaneously optimizes the minimum free energy and codon usage efficiency of the CDS region of an mRNA sequence, this method can take into account the minimum free energy of the full-length mRNA sequence (5'UTR + CDS + 3'UTR) and the codon usage efficiency of the CDS region.
[0115] Figure 2 It is a comparison schematic diagram of the present invention and the prior art.
[0116] Specifically, to solve the problem of optimizing the design of mRNA sequences in the context of 5'UTR and 3'UTR, in this preferred embodiment, a new objective function is first defined. By using the growing method, new codons are added one by one, and the minimum free energy and codon usage efficiency are calculated and scored for ranking. Then, the elimination method is used to exclude sequences containing restriction enzyme sites, bad properties, or poor scores. The clustering search method is used to select representative sequences to increase diversity, and finally, sequences are obtained for recommended synthesis. The sequence design method of this preferred embodiment has the advantages of high sequence innovation degree and fast execution speed, and can be widely used in the task of optimizing the design of mRNA sequences.
[0117] The following will be described in detail:
[0118] 1. Objective Function
[0119] Different from the prior art, in this preferred embodiment, a new objective function (Formula 1) is defined when designing and optimizing mRNA sequences. In this function, the total score of the sequence is the weighted sum of the minimum free energy MFE of the whole sequence whole and the codon usage efficiency CAI of the CDS region CDS . Where K is the weight coefficient. K = 0 means that the objective optimization function only considers the minimum free energy of the full-length sequence; K = +∞ means that the objective optimization function only considers the codon usage efficiency of the CDS region; if K is a positive number, for example, K = 10.0, it means that the objective optimization function needs to balance between the minimum free energy of the whole sequence and the codon usage efficiency of the CDS region. Regarding the calculation method of CAI CDS : First, the CAI value of each codon in the CDS region can be obtained. Secondly, after taking the negative logarithm transformation of this codon CAI value and adding them up, the codon usage efficiency of the CDS region can be obtained (Formula 2).
[0120] Formula 1: Score = MFE whole + K * CAI CDS
[0121] Formula 2:
[0122] 2. Growing Algorithm
[0123] Given a protein sequence, the number of theoretically corresponding codon permutation and combination spaces is extremely large. To avoid exhausting computing power by enumerating all possible combinations, we use the growing method to search for optimized sequences. Figure 3It is a flowchart of this preferred embodiment, and its calculation steps are as follows: (1) At the beginning of the calculation, a given 5'UTR is obtained as the starting input, and then all possible codons corresponding to a residue (amino acid) are added based on this sequence, and the generated new sequences are put into the sequence pool. (2) At this time, judge the number of sequences in the sequence pool: if the number is greater than the search depth D (usually D is set to 100), then proceed to step 3 for calculation; if the number is less than the search depth, return to step 1. (3) Calculate the full-length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS . (4) According to the pre-set target optimization function weight coefficient K value, calculate the score of each sequence and sort them from low to high. (5) According to the set filtering rules or dropout rate, exclude the sequences that do not meet the requirements or discard the sequences with the lowest ranking, and retain the remaining sequences. (6) Among the retained sequences, cluster them according to the search depth. If D = 100, cluster them into 1 hundred classes (clusters), and then select the sequence with the best score from each cluster and put it back into the sequence pool, and eliminate the remaining sequences. And so on, repeat steps 1-6 until all the codons encoding proteins are added. Finally, add the 3'UTR sequence and perform steps 3-4 of the calculation again, and recommend the designed sequence after sorting.
[0124] 3. Dropping
[0125] In the actual synthesis process of mRNA, it is necessary to consider that some sequence fragments (motifs) or sequence properties may increase the synthesis difficulty or cause failure. Therefore, it is the best choice to filter and exclude them during the sequence growth process. For example, the restriction enzyme sites of some DNA plasmids need to be excluded during the design process to avoid incorrect cleavage. Another example is that if the GC content of a certain window of the sequence is too high, it will make the synthesis more difficult. These filtering indicators can all be excluded during the design process. In addition, in order to ensure the convergence of the results, end dropping is also required: that is, in each round of "growth - calculation" cycle, the sequences with poor scores are eliminated so that the calculation results converge to the sequences with better scores.
[0126] 4. Beam Cluster Search
[0127] Since the number of sequences in the sequence pool increases several times (1 to 6 times) in each round of "growth - calculation", in order to avoid an uncontrollable increase in the computational amount, it is necessary to select an appropriate number of sequences to enter the next round of loop. First, we define a search depth D, and D can be equal to 50, 100, or even higher values. The higher the D value, the more sequences can enter the next round, and the greater the computational consumption; vice versa. After multiple steps of growth, the number of sequences in the sequence pool generated in the previous step is fixed at D. After growing each codon, the number of sequences will become D to 6D. Therefore, it is necessary to select D sequences to enter the next round of loop. Beam Search can be used to directly select the top D sequences with the highest scores to enter the next round. However, doing so will make the diversity of the final result worse, that is, the similarity degree among the designed sequences is too high. To avoid this situation, this preferred embodiment adopts a new strategy, that is, clustering the remaining candidate sequences in the sequence pool after elimination into D categories, and then selecting a representative sequence from each category to enter the next round. We name this new method Beam Cluster Search, which can effectively solve the problem of insufficient diversity of the final sequences.
[0128] 5. Recommend
[0129] Limited by the synthesis cost and testing cost, not all of the designed sequences will be synthesized. Therefore, it is very important to recommend some representative sequences. In the last step after adding the 3'UTR sequence, all sequences are scored and sorted, and then clustered. Sequence clustering is performed according to the number of selected sequences R, and then the sequence with the highest score is selected from each category for recommendation. This can effectively ensure the diversity of the recommended sequences.
[0130] Beneficial effects
[0131] This preferred embodiment relates to a method for optimizing messenger ribonucleic acid sequences. In this application, different from the traditional method that only unilaterally considers the properties of the CDS region, this preferred embodiment creatively proposes an mRNA sequence optimization design method in the context of 5'UTR and 3'UTR. This method combines the minimum free energy MFE of the full - length mRNA sequence whole and the codon usage efficiency CAI of the CDS region CDSIntegrated into an objective function, these two objectives can be optimized simultaneously, ensuring the stability and expressivity of the mRNA sequence. In addition, by using the growth method, the computational amount is effectively controlled, avoiding computing power exhaustion. Furthermore, through the elimination method, defective sequences or sequences with poor scores can be avoided, ensuring the quality of the designed sequences. Finally, through beam cluster search and recommendation, the diversity of the designed sequences is improved, enhancing the innovation of the designed sequences.
[0132] Next, a specific embodiment of this application will be introduced.
[0133] The spike protein is a large structural protein on the surface of the coronavirus. It can not only help the virus invade host cells but is also an important immune antigen. The full length of the spike protein of the coronavirus is 1273 amino acids. Both the mRNA vaccines of Moderna and BioNTech / Pfizer encode the spike protein of the coronavirus. After the vaccine is injected into the human body, a large amount of spike protein is produced with the help of human cells, inducing an immune response. After the human body generates immune memory, it can produce a rapid immune response to subsequent possible virus infections, thus avoiding infection or alleviating symptoms.
[0134] From the perspective of mRNA sequence design, the vaccines of Moderna and BioNTech use different 5'UTR and 3'UTR sequences. At the same time, due to the use of different design algorithms by the two companies, the mRNA sequences in the CDS region are also different. Through analysis using the online Clustal Omega program, it is found that the identity of the CDS region sequences of these two vaccines is 90.48%. To verify that the algorithm of the present invention can effectively optimize the minimum free energy (MFE) and codon adaptation index (CAI) of the full-length mRNA sequence, this specific embodiment uses the 5'UTR and 3'UTR sequences of Moderna vaccine mRNA-1273 as constraints to design the mRNA sequence of the spike protein.
[0135] The results show that whether it is the beam search method or the beam cluster search method, compared with Moderna vaccine mRNA-1273, both can improve its minimum free energy and codon adaptation index simultaneously (see Figure 4 and Figure 5 ). This means that the stability and expressivity of the newly designed mRNA sequence can both be increased.
[0136] However, by comparing the identities of the CDS regions of the 10 sequences generated by the beam search method and the beam cluster search method, it can be found that: the identity of the CDS region of the sequences generated by the beam search method is too high (99.8%); on the contrary, without significantly sacrificing the MFE and CAI metrics, the identity of the CDS region of the sequences generated by the beam cluster search method (94.6%) has been greatly improved (as shown in Table 1 below and Figure 6 shown, where, Identity CDS refers to the identity of the CDS region, and Average Identity CDS refers to the average identity of the CDS region). If measured by sequence diversity, the sequence diversity generated by the beam cluster method is more than 46 times that generated by the beam search method.
[0137] Table 1. Design of Spike protein mRNA including 5'UTR and 3'UTR using the beam search (Beam Search) and beam cluster search (Beam Cluster Search) methods (D = 100, R = 10).
[0138]
[0139] Taking the 10 recommended sequences generated by the beam cluster search method (K = 4, D = 100, drop rate = 0.02) as an example. The average identity of the CDS regions of these 10 sequences with the CDS region of the BioNTech vaccine is 87.6%; while the average identity with the CDS region of the Moderna vaccine is 92.1%; the identity between the sequences is 94.6% (see Figure 7 ). In contrast, for the 10 recommended sequences generated by the beam search method (K = 4, D = 100), the average identity of the CDS regions with the CDS region of the BioNTech vaccine is 85.7%; while the average identity with the CDS region of the Moderna vaccine is 88.6%; the identity between the sequences is 99.8% (see Figure 8 ).
[0140] In addition, we also studied the influence of the two parameters of the drop rate and the search depth D in the algorithm on the diversity of the designed sequences. It was found that during the process of increasing the drop rate from 1% to 5%, the diversity of the designed sequences decreased, so choosing a smaller drop rate is beneficial to increasing sequence diversity (see Figure 9 , where, Identity CDS refers to the identity of the CDS region, and Average Identity CDS refers to the average identity of the CDS region). And during the process of increasing the search depth from 50 to 500, the diversity of the designed sequences increased slightly, but not significantly (seeFigure 10 , where Identity CDS refers to the CDS region identity, and Average Identity CDS refers to the average identity of the CDS region).
[0141] In summary, the present preferred embodiment relates to a new method for designing messenger ribonucleic acid (mRNA) sequences. This method can take into account both the minimum free energy of the full-length mRNA sequence (5’UTR + CDS + 3’UTR) and the codon usage efficiency of the CDS region. In this way, its stability can be optimized from the perspective of the overall molecule, and at the same time, the expressibility of the CDS region is optimized. Practice shows that the mRNA sequences designed using this method have better stability and expressibility compared to the mRNA sequences designed using traditional methods. This method can be applied to sequence design and optimization in the fields of mRNA therapeutic drugs, mRNA vaccines, gene therapy drugs, etc.
[0142] Each method implementation manner of the present invention can be implemented in software, hardware, firmware, etc. Whether the present invention is implemented in software, hardware, or firmware, the instruction code can be stored in any type of computer-accessible memory (such as permanent or modifiable, volatile or non-volatile, solid-state or non-solid-state, fixed or replaceable media, etc.). Similarly, the memory can be, for example, Programmable Array Logic (PAL), Random Access Memory (RAM), Programmable Read Only Memory (PROM), Read-Only Memory (ROM), Electrically Erasable Programmable ROM (EEPROM), magnetic disk, optical disk, Digital Versatile Disc (DVD), etc.
[0143] The second embodiment of the present invention relates to a messenger ribonucleic acid sequence design device. Figure 11 It is a schematic structural diagram of the messenger ribonucleic acid sequence design device.
[0144] Specifically, as Figure 11 shown, the messenger ribonucleic acid sequence design device includes:[[]]
[0145] A first sequence generation unit, configured to add, at the end of a given messenger ribonucleic acid sequence, one codon from all possible codons corresponding to the first amino acid of the encoded protein, to generate a set of initial sequences, and place the initial sequences into a sequence pool; wherein, only one codon is added to the end of the sequence, and each codon is added only once;
[0146] A first judgment unit, configured to judge whether the number of sequences in the sequence pool is greater than the search depth D;
[0147] A second sequence generation unit, configured to, when the first judgment unit judges that the number of sequences in the sequence pool is not greater than the search depth D, add, at the end of all sequences in the sequence pool, one codon from all possible codons corresponding to the next amino acid of the encoded protein, so as to generate a set of new sequences; wherein, only one codon is added to the end of the sequence, and each codon is added only once;
[0148] A first calculation unit, configured to, when the first judgment unit judges that the number of sequences in the sequence pool is greater than the search depth D, calculate the full-length minimum free energy MFE of each sequence in the sequence pool whole and the codon adaptation index CAI of the CDS region CDS ;
[0149] A second calculation unit, configured to calculate the score of each sequence in the sequence pool according to the full-length minimum free energy MFE whole output by the first calculation unit CDS and the codon adaptation index CAI of the CDS region, and sort all sequences in the sequence pool according to the scores;
[0150] A filtering unit, configured to filter out sequences that do not meet the requirements in the sequence pool or eliminate sequences with low rankings in the sequence pool according to a preset filtering rule or elimination rate, and retain the remaining sequences;
[0151] A first clustering unit, configured to cluster the retained sequences into D categories according to the search depth D, and retain the sequence with the highest score in each category in the sequence pool, and eliminate the remaining sequences;
[0152] A second judgment unit, configured to judge whether all codons corresponding to all amino acids of the encoded protein have been added;
[0153] The second sequence generation unit is further configured to, when the second determination unit determines that the codons corresponding to all amino acids of the encoded protein have not been added completely, add one codon from all possible codons corresponding to the next amino acid of the encoded protein to the end of each sequence in the sequence pool, so as to generate a group of new sequences; wherein, only one codon is added to the end of the sequence, and each codon is added only once;
[0154] The 3' UTR sequence addition unit is configured to, when the second determination unit determines that the codons corresponding to all amino acids of the encoded protein have been added completely, add 3' UTR sequences to the end of each sequence in the sequence pool;
[0155] The third calculation unit is configured to calculate the full-length minimum free energy MFE of each sequence in the sequence pool after adding the 3' UTR sequence whole and the codon adaptation index CAI of the CDS region CDS ;
[0156] The fourth calculation unit is configured to calculate the score of each sequence in the sequence pool after adding the 3' UTR sequence according to the full-length minimum free energy MFE whole output by the third calculation unit and the codon adaptation index CAI of the CDS region CDS , and sort all the sequences in the sequence pool after adding the 3' UTR sequence according to the scores;
[0157] The second clustering unit is configured to select the number R of sequences, cluster all the sequences in the sequence pool after adding the 3' UTR sequence into R categories, and recommend the sequence with the highest score in each category as the designed messenger ribonucleic acid sequence.
[0158] The embodiment of the present application provides a method for designing a messenger ribonucleic acid (mRNA) sequence in an environment with 5' UTR and 3' UTR, which can simultaneously optimize the full-length minimum free energy of an mRNA sequence and the codon adaptation index of the CDS region, ensuring the stability and expressibility of the mRNA sequence.
[0159] This embodiment is a device embodiment corresponding to the first embodiment, and this embodiment can be implemented in cooperation with the first embodiment. The relevant technical details mentioned in the first embodiment are still valid in this embodiment. To avoid repetition, they are not elaborated here. Correspondingly, the relevant technical details mentioned in this embodiment can also be applied in the first embodiment.
[0160] It should be noted that each unit mentioned in the embodiments of each device of the present invention is a logical unit. Physically, a logical unit can be a physical unit, a part of a physical unit, or can be implemented as a combination of multiple physical units. The physical implementation manner of these logical units themselves is not the most important. The combination of the functions implemented by these logical units is the key to solving the technical problems proposed by the present invention. In addition, in order to highlight the innovative part of the present invention, the above embodiments of each device of the present invention do not introduce units that are not closely related to solving the technical problems proposed by the present invention. This does not mean that there are no other units in the above device embodiments.
[0161] It should be noted that those skilled in the art should understand that the implementation functions of each unit shown in the above embodiments of each device can be understood with reference to the relevant descriptions of the corresponding methods described above. The functions of each unit shown in the above embodiments of each device can be implemented by a program (executable instruction) running on a processor or by specific logic circuits. If the above devices of the embodiments of this specification are implemented in the form of software function modules and sold or used as independent products, they can also be stored in a computer-readable storage medium. Based on such an understanding, the technical solutions of the embodiments of this specification, in essence, or the part that contributes to the prior art, can be embodied in the form of a software product. This computer software product is stored in a storage medium and includes several instructions for causing a computer device (which can be a personal computer, a server, or a network device, etc.) to execute all or part of the methods described in the embodiments of this specification. The aforementioned storage medium includes: various media that can store program codes such as USB flash drives, mobile hard disks, read-only memories (ROMs), magnetic disks, or optical discs. In this way, the embodiments of this specification are not limited to any specific combination of hardware and software.
[0162] Accordingly, an embodiment of the present application further provides a messenger ribonucleic acid sequence design computing device, which includes a memory for storing computer-executable instructions, and a processor; the processor is configured to implement the steps in the above method embodiments when executing the computer-executable instructions in the memory. Among them, the processor may be a central processing unit (Central Processing Unit, abbreviated as "CPU"), or other general-purpose processors, digital signal processors (Digital Signal Processor, abbreviated as "DSP"), application specific integrated circuits (Application Specific Integrated Circuit, abbreviated as "ASIC"), etc. The aforementioned memory may be a read-only memory (read-only memory, abbreviated as "ROM"), random access memory (random access memory, abbreviated as "RAM"), flash memory (Flash), hard disk or solid state drive, etc. The steps of the methods disclosed in the embodiments of the present invention can be directly implemented by a hardware processor, or implemented by a combination of hardware and software modules in the processor.
[0163] In addition, an embodiment of the present application further provides a computer storage medium, in which computer-executable instructions are stored, and the computer-executable instructions implement the method embodiments of the present application when executed by a processor.
[0164] It should be noted that in the application documents of this patent, relational terms such as first and second are only used to distinguish one entity or operation from another entity or operation, and do not necessarily require or imply any such actual relationship or order between these entities or operations. Moreover, the term "comprising", "including" or any other variation thereof is intended to cover non-exclusive inclusion, so that a process, method, article or device including a series of elements includes not only those elements, but also other elements not expressly listed, or elements inherent to such process, method, article or device. Without further limitation, an element defined by the statement "including one" does not exclude the existence of another identical element in the process, method, article or device including the element. In the application documents of this patent, if it is mentioned that an act is performed according to a certain element, it means at least performing the act according to the element, including two cases: performing the act only according to the element, and performing the act according to the element and other elements. Expressions such as multiple, multiple times, multiple types, etc. include 2, 2 times, 2 types, as well as more than 2, more than 2 times, more than 2 types.
[0165] All documents mentioned in this application are considered to be included in the disclosure of this application as a whole so that they can be used as a basis for modification if necessary. In addition, it should be understood that after reading the above disclosure of this application, those skilled in the art can make various changes or modifications to this application, and these equivalent forms also fall within the scope of protection required by this application.
Claims
1. A method for designing a messenger ribonucleic acid sequence, characterized in that, Including the following steps: Step S1: Add one of all possible codons corresponding to the first amino acid of the encoded protein to the end of the given messenger ribonucleic acid sequence respectively to generate a set of initial sequences, and put the initial sequences into a sequence pool; wherein, only one codon is added to the end of the sequence, and each codon is only added once; Step S2: Determine whether the number of sequences in the sequence pool is greater than the search depth D; If not, then execute Step S3: Add one of all possible codons corresponding to the next amino acid of the encoded protein to the end of all sequences in the sequence pool respectively, so as to generate a set of new sequences, and return to execute Step S2; wherein, only one codon is added to the end of the sequence, and each codon is only added once; If so, perform step S4: calculate the full-length minimum free energy MFE of each sequence in the sequence pool whole and the codon adaptation index CAI of the CDS region CDS ; Step S5: According to the full-length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS , calculate the score of each sequence in the sequence pool, and sort all the sequences in the sequence pool according to the score; Step S6: Filter out the sequences that do not meet the requirements in the sequence pool or eliminate the sequences with lower rankings in the sequence pool according to the preset filtering rules or elimination rate, and the remaining sequences are retained; Step S7: Cluster the retained sequences into D categories according to the search depth D, and retain the sequence with the highest score in each category in the sequence pool, and eliminate the remaining sequences; Step S8: Determine whether all codons corresponding to all amino acids of the encoded protein have been added; If not, then return to execute Step S3; if so, then execute Step S9: Add the 3'UTR sequence to the end of all sequences in the sequence pool respectively; Step S10: Calculate the full-length minimum free energy MFE of each sequence in the sequence pool after adding the 3'UTR sequence whole and the codon adaptation index CAI of the CDS region CDS ; Step S11: Calculate the score of each sequence in the sequence pool after adding the 3'UTR sequence according to the full-length minimum free energy MFE whole and the codon adaptation index CAI of the CDS region CDS , and sort all the sequences in the sequence pool after adding the 3'UTR sequence according to the score; Step S12: Cluster all sequences after adding the 3'UTR sequence in the sequence pool into R categories according to the number R of selected sequences, and recommend the sequence with the highest score in each category as the designed messenger ribonucleic acid sequence.
2. The method according to claim 1, characterized in that, The score of the said sequence is the full-length minimum free energy MFE of the whole sequence whole and the codon usage efficiency CAI of the said CDS region CDS is a weighted sum, and the formula for the score of the said sequence is expressed as: Score=MFE whole +K*CAI CDS Wherein, Score is the score of the sequence, and K is the weight coefficient.
3. The method according to claim 2, wherein The formula for the score of the sequence is further expressed as: Wherein, Score is the score of the sequence, K is the weight coefficient, m is the number of codons in the CDS region, and i represents the i-th codon in the CDS region.
4. The method according to claim 2 or 3, characterized in that, The value range of K is [0, +∞].
5. The method according to claim 1, wherein In Step S1, the given messenger ribonucleic acid sequence is the 5'UTR sequence.
6. The method according to claim 1, characterized in that, The number of all possible codons corresponding to one amino acid of the encoded protein is 1-6.
7. The method according to claim 1, characterized in that The search depth D is any one of 50, 100, 200, and 500.
8. A messenger ribonucleic acid sequence design device, characterized in that, Including: The first sequence generation unit is used to add one of all possible codons corresponding to the first amino acid of the encoded protein to the end of the given messenger ribonucleic acid sequence respectively to generate a set of initial sequences, and put the initial sequences into a sequence pool; wherein, only one codon is added to the end of the sequence, and each codon is only added once; The first judgment unit is used to judge whether the number of sequences in the sequence pool is greater than the search depth D; A second sequence generation unit, configured to, when the first determination unit determines that the number of sequences in the sequence pool is not greater than the search depth D, add, to the end of each sequence in the sequence pool, one of all possible codons corresponding to the next amino acid encoding the protein, so as to generate a set of new sequences; wherein, only one codon is added to the end of the sequence, and each codon is added only once; A first calculation unit, configured to calculate the full-length minimum free energy MFE of each sequence in the sequence pool when the first determination unit determines that the number of sequences in the sequence pool is greater than the search depth D whole and the codon usage efficiency CAI of the CDS region CDS ; A second computing unit, configured to calculate the score of each sequence in the sequence pool according to the full-length minimum free energy MFE output by the first computing unit whole and the codon adaptation index CAI of the CDS region CDS and sort all the sequences in the sequence pool according to the score; A filtering unit, configured to filter out sequences that do not meet the requirements in the sequence pool or eliminate sequences with a lower ranking in the sequence pool according to a preset filtering rule or elimination rate, and the remaining sequences are retained; A first clustering unit, configured to cluster the retained sequences into D categories according to the search depth D, and retain the sequence with the highest score in each category in the sequence pool, and eliminate the remaining sequences; A second determination unit, configured to determine whether all codons corresponding to all amino acids of the encoded protein have been added; The second sequence generation unit is further configured to, when the second determination unit determines that not all codons corresponding to all amino acids of the encoded protein have been added, add, to the end of each sequence in the sequence pool, one of all possible codons corresponding to the next amino acid encoding the protein, so as to generate a set of new sequences; wherein, only one codon is added to the end of the sequence, and each codon is added only once; A 3'UTR sequence addition unit, configured to, when the second determination unit determines that all codons corresponding to all amino acids of the encoded protein have been added, add 3'UTR sequences to the end of each sequence in the sequence pool; A third calculation unit, configured to calculate the full-length minimum free energy MFE of each sequence in the sequence pool after adding the 3'UTR sequence whole and the codon adaptation index CAI of the CDS region CDS ; A fourth calculation unit, configured to calculate the score of each sequence in the sequence pool after adding the 3'UTR sequence according to the full-length minimum free energy MFE output by the third calculation unit whole and the codon adaptation index CAI of the CDS region CDS , and sort all the sequences in the sequence pool after adding the 3'UTR sequence according to the score; A second clustering unit, configured to select the number R of sequences, cluster all sequences in the sequence pool after adding 3'UTR sequences into R categories, and recommend the sequence with the highest score in each category as the designed messenger ribonucleic acid sequence.
9. A messenger ribonucleic acid sequence design computing device, characterized in that, Comprising: A memory, configured to store computer-executable instructions; And, A processor, configured to implement the steps in the method according to any one of claims 1 to 7 when executing the computer-executable instructions.
10. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores computer-executable instructions, and when the computer-executable instructions are executed by the processor, the steps in the method according to any one of claims 1 to 7 are implemented.
Citation Information
Patent Citations
Neoantigen identification using hotspots
CN111465989A
Method, device and equipment for optimizing 5'untranslated region sequence of messenger ribonucleic acid
CN116168764A