Method and device for designing messenger ribonucleic acid sequence, computing equipment and storage medium

Abstract

The invention relates to the technical field of biological pharmacy, and discloses a messenger ribonucleic acid sequence design method and device, computing equipment and a storage medium. According to the method, on the basis that a given 5 'UTR sequence serves as initial input, new codons are added one by one through a growth method, the overall-length minimum free energy and codon use efficiency of the whole sequence are calculated, scored and sorted, an elimination method is used for excluding sequences containing restriction enzyme cutting sites / poor properties / poorer scores, and finally, the sequences containing the restriction enzyme cutting sites / poor properties / poorer scores are obtained. And selecting representative sequences by using a cluster clustering search method to increase diversity, and finally obtaining the sequences for recommendation synthesis. The invention provides an mRNA sequence optimization design method in an environment with 5 'UTR and 3' UTR, the full-length minimum free energy of an mRNA sequence and the codon use efficiency of a CDS region are considered, the two targets can be optimized at the same time, the stability can be optimized from the perspective of overall molecules, the expressibility of the CDS region is optimized at the same time, and the method has the advantages of being simple in structure, convenient to operate and the like. And the stability and expressibility of the mRNA sequence are ensured.

Description

Technical Field

[0001] The present invention relates to the field of biopharmaceutical technology, and in particular to a messenger RNA sequence design method, device, computing equipment and storage medium. Background Art

[0002] Messenger RNA (mRNA) is a single-stranded molecule that is responsible for transferring genetic information from deoxyribonucleic acid (DNA) to the ribosome, which then decodes the genetic information and synthesizes proteins. In recent years, the technical means of synthesizing mRNA through in vitro transcription (IVT) technology has gradually matured. Biopharmaceutical companies use mRNA technology platforms and targeted, efficient, and low-immunogenic delivery vectors to deliver mRNA into the human body, and human cells use synthetic mRNA to express specific proteins. Because synthetic mRNA is similar in structure to natural mRNA, it enables patients to produce therapeutic proteins in their own bodies, thereby reducing the various troubles caused by the complex manufacturing process of recombinant proteins. Based on the diversity of mRNA therapies, new drug developers use synthetic mRNA as a disease treatment tool and have made progress in the development of cancer immunotherapy, gene therapy, and control of infectious diseases. mRNA technology has obvious efficiency and protection advantages in vaccine development. Moderna and BioNTech both developed and marketed mRNA vaccines within a year, making history.

[0003] When designing mRNA sequences, biopharmaceutical companies often do not directly use wild-type mRNA sequences, but optimize the design of mRNA sequences based on activity and drugability. Traditional mRNA sequence design methods usually only consider codon optimization, that is, using codons with higher frequency of use in human cells to optimize the CDS sequence encoding proteins, so that the synthetic mRNA sequence has higher expressibility than the wild-type mRNA sequence. However, this design method that only considers the frequency of codon usage lacks consideration of physicochemical properties, making the designed synthetic mRNA molecules unstable and easily degraded in vitro. Therefore, low temperature (-20 / -80°C) storage conditions are required, which increases storage and transportation costs.

[0004] In response to the above problems, Huang Liang and others from Baidu Research Institute in the United States proposed a new algorithm for mRNA sequence design in 2020: Linear Design. In this algorithm, the optimization of the CDS region of synthetic mRNA is set as a multi-objective optimization problem: that is, the minimum free energy (MFE) and codon usage efficiency (CAI) of the CDS region are optimized simultaneously. By defining a scoring function S = MFE CDS +λ|p|(-logCAI CDS ), Linear Design simultaneously optimizes the minimum free energy and codon usage efficiency of the CDS region of an mRNA sequence to obtain a synthetic mRNA molecule with good expression and stability. However, in addition to the CDS region encoding proteins, the synthetic mRNA molecules used in the pharmaceutical industry also include non-coding regions such as 5'UTR and 3'UTR, which also have an important impact on the expression and stability of mRNA molecules.

[0005] Therefore, there is an urgent need for a messenger RNA sequence design method that can simultaneously optimize the full-length minimum free energy of an mRNA sequence and the codon usage efficiency of the CDS region. Summary of the invention

[0006] The purpose of the present invention is to provide a messenger RNA sequence design method, device, computing equipment and storage medium, which can simultaneously optimize the full-length minimum free energy of an mRNA sequence and the codon usage efficiency of the CDS region, thereby ensuring the stability and expressibility of the mRNA sequence.

[0007] In order to solve the above technical problems, the embodiments of the present invention disclose a method for designing a messenger RNA sequence, comprising the following steps:

[0008] Step S1: adding one codon from all possible codons corresponding to the first amino acid encoding the protein to the end of a given messenger RNA sequence, generating a set of initial sequences, and placing 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;

[0009] Step S2: determining 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 encoding the protein to the end of all sequences in the sequence pool, thereby generating 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 added only once;

[0011] If yes, then execute step S4: calculate the full-length minimum free energy MFE of each sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ;

[0012] Step S5: According to the full-length minimum free energy MFE whole and the codon usage efficiency CAI of the CDS region CDS , calculating a score for each sequence in the sequence pool, and sorting all sequences in the sequence pool according to the score;

[0013] Step S6: According to a preset filtering rule or elimination rate, the sequences in the sequence pool that do not meet the requirements are filtered out or the sequences with the lowest ranking in the sequence pool are eliminated, and the remaining sequences are retained;

[0014] Step S7: clustering the retained sequences into D categories according to the search depth D, retaining the sequence with the highest score in each category in the sequence pool, and eliminating the remaining sequences;

[0015] Step S8: determining whether all the codons corresponding to the amino acids of the encoded protein have been added;

[0016] If not, return to step S3; if yes, execute step S9: add 3'UTR sequences to the end of all sequences in the sequence pool;

[0017] Step S10: Calculate the full-length minimum free energy MFE of each sequence after adding the 3'UTR sequence to the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ;

[0018] Step S11: According to the full-length minimum free energy MFE whole and the codon usage efficiency CAI of the CDS region CDS , calculating the score of each sequence after the 3'UTR sequence is added to the sequence pool, and sorting all the sequences after the 3'UTR sequence is added to the sequence pool according to the score;

[0019] Step S12: according to the number R of recommended sequences, all sequences in the sequence pool after adding the 3'UTR sequence are clustered into R categories, and the sequence with the highest score in each category is recommended as the designed messenger RNA sequence.

[0020] The embodiment of the present invention also discloses a device for designing RNA sequences, comprising:

[0021] The first sequence generation unit is used to add one codon from all possible codons corresponding to the first amino acid encoding the protein at the end of a given messenger RNA sequence, 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 judging unit, configured to judge whether the number of sequences in the sequence pool is greater than a search depth D;

[0023] A second sequence generating unit is used for adding one codon from all possible codons corresponding to the next amino acid encoding the protein to the end of each sequence in the sequence pool when the first judging unit judges that the number of sequences in the sequence pool is not greater than the search depth D, thereby generating a group of new sequences; wherein only one codon is added to the end of the sequence, and each codon is added only once;

[0024] The first calculation unit is used to calculate the full-length minimum free energy MFE of each sequence in the sequence pool when the first judgment unit determines that the number of sequences in the sequence pool is greater than the search depth D. whole and codon usage efficiency CAI in the CDS region CDS ;

[0025] A second calculation unit is used to calculate the full-length minimum free energy MFE output by the first calculation unit. whole and the codon usage efficiency CAI of the CDS region CDS , calculating a score for each sequence in the sequence pool, and sorting all sequences in the sequence pool according to the score;

[0026] A filtering unit, used 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;

[0027] A first clustering unit is used 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 is used to judge whether the codons corresponding to all the amino acids of the encoded protein have been added;

[0029] The second sequence generating unit is further configured to, when the second judging unit judges that the codons corresponding to all the amino acids of the encoded protein have not been added, add one codon from all the possible codons corresponding to the next amino acid of the encoded protein to the end of all the sequences in the sequence pool, thereby generating a set of new sequences; wherein only one codon is added to the end of the sequence, and each codon is added only once;

[0030] A 3'UTR sequence adding unit, used for adding 3'UTR sequences to the ends of all sequences in the sequence pool respectively when the second judgment unit judges that the codons corresponding to all the amino acids of the encoded protein have been added;

[0031] The third calculation unit is used to calculate the full-length minimum free energy MFE of each sequence after adding the 3'UTR sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ;

[0032] The fourth calculation unit is used to calculate the full-length minimum free energy MFE output by the third calculation unit. whole and the codon usage efficiency CAI of the CDS region CDS , calculating the score of each sequence after the 3'UTR sequence is added to the sequence pool, and sorting all the sequences after the 3'UTR sequence is added to the sequence pool according to the score;

[0033] The second clustering unit is used to recommend the number R of sequences, cluster all the sequences after the 3'UTR sequence is added to the sequence pool into R categories, and recommend the sequence with the highest score in each category as the designed messenger RNA sequence.

[0034] The embodiment of the present invention further discloses a messenger RNA sequence design computing device, comprising:

[0035] a memory for storing computer executable instructions; and,

[0036] A processor is used to implement the steps in the above method when executing the computer executable instructions.

[0037] The embodiment of the invention further discloses a computer-readable storage medium, wherein the computer-readable storage medium stores computer-executable instructions, and the computer-executable instructions implement the steps in the above method when executed by a processor.

[0038] Compared with the prior art, the main differences and effects of the embodiments of the present invention are:

[0039] The full-length minimum free energy of an mRNA sequence and the codon usage efficiency of the CDS region can be optimized simultaneously, ensuring the stability and expressibility of the mRNA sequence.

[0040] Furthermore, by using the growth method, the amount of calculation is effectively controlled and computing power exhaustion is avoided.

[0041] Furthermore, through the elimination method, defective sequences or sequences with poor scores can be avoided to ensure the quality of the designed sequence.

[0042] Furthermore, through clustering search and recommendation, the diversity of design sequences is improved and the innovation of design sequences is enhanced.

[0043] A large number of technical features are recorded in the specification of this application, which are distributed in various technical solutions. If all possible combinations of technical features of this application (i.e., technical solutions) are to be listed, the specification will be too long. In order to avoid this problem, the various technical features disclosed in the above-mentioned invention content of this application, the various technical features disclosed in the various embodiments and examples below, and the various technical features disclosed in the accompanying drawings can be freely combined with each other to form various new technical solutions (these technical solutions are all deemed to have been recorded in this specification), unless the combination of such technical features is technically infeasible. For example, in one example, feature A+B+C is disclosed, and in another example, feature A+B+D+E is disclosed, and features C and D are equivalent technical means that play the same role. Technically, only one can be used, and it is impossible to use them at the same time. Feature E can be combined with feature C technically. Then, the solution of A+B+C+D should not be deemed to have been recorded because it is technically infeasible, and the solution of A+B+C+E should be deemed to have been recorded. BRIEF DESCRIPTION OF THE DRAWINGS

[0044] Figure 1 is a schematic flow chart of a messenger RNA sequence design method in the first embodiment of the present invention;

[0045] Figure 2 is a schematic diagram comparing the first embodiment of the present invention with the prior art;

[0046] Figure 3 It is a schematic flow chart of a method for designing a messenger RNA sequence in a preferred embodiment of the present invention;

[0047] Figure 4 It is to optimize the objective function MFE whole +K*CAI CDSBelow is a schematic diagram of using the beam search method to design mRNA for Spike protein under specific UTR constraints;

[0048] Figure 5 It is to optimize the objective function MFE whole +K*CAI CDS Below is a schematic diagram of using the Beam Cluster Search method to design mRNA for Spike protein under specific UTR constraints;

[0049] Figure 6 It is a schematic diagram of using beam search and beam cluster search methods to design the comparison of the CDS region identity of the Spike protein mRNA sequence;

[0050] Figure 7 This is a schematic diagram of the comparison of the similarity between the Spike protein mRNA sequence (K=4, D=100, drop rate=0.02) and the CDS region of the BioNTech and Moderna marketed vaccines designed by the Beam Cluster Search method;

[0051] Figure 8 This is a schematic diagram of the comparison of the identity between the Spike protein mRNA sequence (K=4, D=100) and the CDS region of the BioNTech and Moderna marketed vaccines designed using the beam search method;

[0052] Fig. 9 It is a schematic diagram of the effect of the drop rate (DR) on the diversity of the designed sequence;

[0053] Fig.10 It is a schematic diagram of the effect of search depth (Depth, D) on the diversity of designed sequences;

[0054] Fig.11 It is a schematic diagram of the structure of a messenger RNA sequence design device in the second embodiment of the present invention. DETAILED DESCRIPTION

[0055] In the following description, many technical details are provided to help readers better understand the present application. However, those skilled in the art can understand that even without these technical details and various changes and modifications based on the following embodiments, the technical solutions claimed in the claims of the present application can be implemented.

[0056] Description of some concepts:

[0057] mRNA: messenger ribonucleic acid

[0058] MFE (minimum free energy): minimum free energy

[0059] CAI (codon adaption index): codon usage efficiency

[0060] diversity: diversity

[0061] In order to make the objectives, technical solutions and advantages of the present invention more clear, 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 RNA sequence. Figure 1 It is a schematic diagram of the process of the messenger RNA sequence design method.

[0063] For a given protein sequence, the number of codon permutations and combinations corresponding to it is theoretically very large. To avoid exhausting computing power by enumerating all possible operations, we use the growing algorithm to search for optimized sequences.

[0064] Specifically, if Figure 1 As shown, the messenger RNA sequence design method comprises the following steps:

[0065] In step S1, one codon from all possible codons corresponding to the first amino acid encoding the protein is added to the end of a given messenger RNA sequence to generate a set of initial sequences, and the initial sequences are placed in a sequence pool; wherein only one codon is added to the end of the sequence, and each codon is added only once.

[0066] In this embodiment, preferably, the given messenger RNA sequence is a 5'UTR sequence.

[0067] That is, in step S1, the starting sequence is the 5'UTR sequence.

[0068] The 5'UTR sequence is used as the starting input, and then based on the 5'UTR sequence, all possible codons corresponding to the first amino acid of the encoded protein are added to the end of the 5'UTR sequence to generate a set of initial sequences, and the initial sequences are placed in the sequence pool.

[0069] It should be noted that only one codon is added to the end of the 5'UTR sequence each time, and each of all possible codons corresponding to the first amino acid of the encoded protein is added only 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] Accordingly, in step S1, one codon is added to the end of the 5'UTR start sequence each time. After all possible codons corresponding to the first amino acid of the encoded protein are added, the 5'UTR start sequence grows into 1-6 initial sequences.

[0072] Then, the process proceeds to step S2: determining whether the number of sequences in the sequence pool is greater than the search depth D.

[0073] If not, execute step S3: add one codon from all possible codons corresponding to the next amino acid encoding the protein to the end of all sequences in the sequence pool, thereby generating 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 added only once.

[0074] If yes, then execute step S4: calculate the full-length minimum free energy MFE of each sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS .

[0075] Then, step S5 is entered: according to the full-length minimum free energy MFE whole and the codon usage efficiency 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.

[0076] In this embodiment, preferably, the score of the sequence is the minimum free energy MFE of the entire sequence. whole Codon usage efficiency 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] Further, 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 step S5, the score of each sequence is calculated according to the preset target optimization weight coefficient K value, and the scores are sorted 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 entire 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 must achieve a balance between the minimum free energy of the entire length of the sequence and the codon usage efficiency of the CDS region.

[0085] Then, the process proceeds to step S6: according to a preset filtering rule or a drop rate, the sequences in the sequence pool that do not meet the requirements are filtered out or the sequences with the lowest ranking in the sequence pool are eliminated, and the remaining sequences are retained.

[0086] In step S6, the elimination method is used to exclude sequences containing restriction sites / bad properties / or poor scores.

[0087] Then, the process proceeds to step S7: clustering the retained sequences into D categories according to the search depth D, retaining the sequence with the highest score in each category in the sequence pool, and eliminating the remaining sequences.

[0088] In this implementation, preferably, the search depth D may be any one of 50, 100, 200 and 500.

[0089] In this implementation, the search depth D is taken as 100 for illustration.

[0090] Of course, this is only a preferred implementation of the present application, and D can be equal to any one of 50, 100, 200 and 500, or even a higher value, and is not limited thereto.

[0091] The higher the D value, the more sequences will enter the next round, the better the diversity of the designed sequences will be, and of course, the greater the computational cost will be.

[0092] In step S7, clustering is performed on the retained sequences according to the search depth. If D=100, the sequences are clustered into 100 clusters. Then, the sequences with the highest scores are selected from each cluster and retained in the sequence pool, and the remaining sequences are eliminated.

[0093] Then, the process proceeds to step S8: determining whether the codons corresponding to all the amino acids of the encoded protein have been added.

[0094] If not, return to step S3.

[0095] And so on, steps S3-S8 are repeated until all protein encoding codons are added.

[0096] After multiple rounds of growth, the number of sequences generated in the previous round that remain in the sequence pool is fixed at D. After each codon is grown, the number of sequences in the sequence pool will become D-6D, so D sequences must be selected to enter the next round of cycles.

[0097] If yes, then execute step S9: add 3'UTR sequences to the end of all sequences in the sequence pool.

[0098] Finally, the 3'UTR sequence is added to grow into a full-length mRNA sequence (5'UTR+CDS+3'UTR). In addition to the CDS region encoding the protein, the full-length mRNA sequence also includes non-coding regions such as 5'UTR and 3'UTR, which also have an important impact on the expression and stability of the mRNA molecule.

[0099] Then, step S10 is entered: calculating the full-length minimum free energy MFE of each sequence after adding the 3'UTR sequence to the sequence pool whole and codon usage efficiency CAI in the CDS region CDS .

[0100] Then, step S11 is entered: according to the full-length minimum free energy MFE whole and the codon usage efficiency CAI of the CDS region CDS , calculating the score of each sequence after adding the 3'UTR sequence to the sequence pool, and sorting all the sequences after adding the 3'UTR sequence to the sequence pool according to the score.

[0101] Similarly, preferably, in step S11, the score of the sequence is the minimum free energy MFE of the entire sequence. whole The codon usage efficiency of the CDS region CDS The formula for the score of the sequence is expressed as:

[0102] Score=MFE whole +K*CAI CDS

[0103] Wherein, 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 step S11, the score of each sequence is calculated according to the preset target optimization weight coefficient K value, and the scores are sorted 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 entire 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 must achieve a balance between the minimum free energy of the entire length of the sequence and the codon usage efficiency of the CDS region.

[0110] Then, step S12 is entered: according to the number R of recommended sequences, all sequences after adding 3'UTR sequences in the sequence pool are clustered into R categories, and the sequence with the highest score in each category is recommended as the designed messenger RNA sequence.

[0111] This process ends thereafter.

[0112] In summary, the embodiments of the present application propose a method for designing messenger RNA (mRNA) sequences in the presence of 5'UTR and 3'UTR, which can simultaneously optimize the full-length minimum free energy of an mRNA sequence and the codon usage efficiency of the CDS region, thereby 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 the preferred embodiment are mainly for ease of understanding and are not intended to limit the scope of protection of this application.

[0114] The technical solution of this preferred embodiment proposes a new mRNA sequence design optimization method. Unlike the prior art, which only 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 and codon usage efficiency of the CDS region of the full-length mRNA sequence (5'UTR+CDS+3'UTR).

[0115] Figure 2 It is a comparative schematic diagram of the present invention and the prior art.

[0116] Specifically, in order to solve the problem of optimizing the design of mRNA sequences in the presence of 5'UTR and 3'UTR, in the preferred embodiment, a new objective function is first defined, and new codons are added one by one by using the growth method and the minimum free energy and codon usage efficiency are calculated and scored, and then the elimination method is used to exclude sequences containing restriction sites / bad properties / or poor scores, and then the clustering search method is used to select representative sequences to increase diversity, and finally the sequence is obtained for recommended synthesis. The sequence design method of the preferred embodiment has the advantages of high sequence innovation and fast execution speed, and can be widely used in mRNA sequence optimization design tasks.

[0117] The following will explain in detail:

[0118] 1. Objective Function

[0119] Different from the prior art, this preferred embodiment defines a new objective function (Formula 1) when designing and optimizing mRNA sequences. In this function, the total score of the sequence is the minimum free energy MFE of the entire sequence. whole Codon usage efficiency CAI CDS The weighted sum of K is the weight coefficient. K = 0 means that the target optimization function only considers the minimum free energy of the full-length sequence; K = +oe means that the target optimization function only considers the codon usage efficiency of the CDS region; if K is a positive number, such as K = 10.0, it means that the target optimization function must achieve a balance between the minimum free energy of the overall sequence and the codon usage efficiency of the CDS region. About CAI CDS Calculation method: First, the CAI value of each codon in the CDS region can be obtained. Secondly, the codon CAI value can be negatively logarithmically transformed and then added to obtain the codon usage efficiency of the CDS region (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 codon permutations and combinations corresponding to it is theoretically very large. To avoid exhausting computing power by enumerating all possible operations, we use the growth method to search for optimized sequences. Figure 3This is a flow chart of the preferred embodiment, and the 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 sequence is placed in the sequence pool. (2) At this time, the number of sequences in the sequence pool is judged: if the number is greater than the search depth D (usually D is set to 100), the third step is performed; if the number is less than the search depth, return to step 1. (3) Calculate the full-length minimum free energy MFE of each mRNA sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS . (4) According to the preset 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 elimination rate (drop rate), eliminate the sequences that do not meet the requirements or discard the sequences with the lowest ranking, and retain the remaining sequences. (6) In the retained sequences, cluster according to the search depth. If D = 100, cluster into 100 classes (cluster), and then select the sequence with the best score from each cluster and put it back into the sequence pool, and the remaining sequences are eliminated. And so on, repeat steps 1-6 until all protein encoding codons have been added. Finally, add the 3'UTR sequence, perform steps 3-4 again, and recommend the output design sequence after sorting.

[0124] 3. Dropping

[0125] In the actual synthesis process of mRNA, it is necessary to consider that certain sequence fragments (motifs) or sequence properties may increase the difficulty of synthesis or cause failure. Therefore, filtering and excluding them during the sequence growth process is the best choice. For example, the restriction sites of certain DNA plasmids need to be excluded during the design process to avoid incorrect shearing. For example, 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 be excluded during the design process. In addition, in order to ensure the convergence of the results, tail elimination 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 each round of "growth-calculation" will increase the number of sequences in the sequence pool by several times (1 to 6 times), in order to avoid an uncontrollable increase in the amount of calculation, an appropriate number of sequences must be selected to enter the next round of cycles. First, we define a search depth D, which can be equal to 50, 100, or even higher. 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 generated in the previous step that enter the sequence pool is fixed to D. After each codon is grown, the number of sequences will become D to 6D, so D sequences must be selected to enter the next round of cycles. Beam Search can be used to directly select the top D sequences in the score ranking to enter the next round. However, doing so will make the diversity of the final result worse, that is, the degree of identity between the designed sequences is too high. To avoid this situation, this preferred embodiment adopts a new strategy, that is, clustering the candidate sequences remaining after elimination in the sequence pool into D categories, and then selecting a representative sequence from each category to enter the next round. We named this new method Beam Cluster Search, which can effectively solve the problem of insufficient final sequence diversity.

[0128] 5. Recommend

[0129] Due to the limitation of synthesis and testing costs, not all designed sequences will be synthesized, so 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. Sequences are clustered according to the number R of recommended sequences, and then a sequence with the highest score is selected from each category for recommendation. This can effectively ensure the diversity of recommended sequences.

[0130] Beneficial Effects

[0131] This preferred embodiment relates to a method for optimizing messenger RNA sequences. In this application, unlike the traditional method that only considers the properties of the CDS region unilaterally, this preferred embodiment creatively proposes a method for optimizing the design of mRNA sequences in the presence of 5'UTR and 3'UTR. This method uses the minimum free energy MFE of the entire length of the mRNA sequence as the whole and codon usage efficiency CAI in the CDS region CDSBy integrating them into one 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 amount of calculation is effectively controlled and the exhaustion of computing power is avoided. Secondly, through the elimination method, defective sequences or sequences with poor scores can be avoided to ensure the quality of the designed sequence. Finally, through clustering search and recommendation, the diversity of the designed sequence is improved and the innovation of the designed sequence is improved.

[0132] A specific embodiment of the present application is introduced below.

[0133] The spike protein is a huge structural protein located on the surface of the coronavirus. It not only helps the virus invade host cells, but is also an important immune antigen. The spike protein of the coronavirus is 1273 amino acids in length. Both Moderna and BioNTech / Pfizer's mRNA vaccines 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 to induce an immune response. After the human body produces immune memory, it can produce a rapid immune response to subsequent possible viral infections, thereby avoiding infection or alleviating symptoms.

[0134] From the perspective of mRNA sequence design, Moderna and BioNTech's vaccines use different 5'UTR and 3'UTR sequences. At the same time, because the two companies used different design algorithms, the mRNA sequences in the CDS region were also different. By using the online Clustal Omega program for analysis, it was found that the sequence identity (Identity) of the CDS region of the two vaccines was 90.48%. In order to verify that the algorithm of the present invention can effectively optimize the minimum free energy (MFE) and codon usage efficiency (CAI) of the full-length mRNA sequence, this specific embodiment uses the 5'UTR and 3'UTR sequences of the Moderna vaccine mRNA-1273 as constraints to design the mRNA sequence for the Spike protein.

[0135] The results show that both the beam search method and the beam cluster search method can improve the minimum free energy and codon usage efficiency of the Moderna vaccine mRNA-1273 (see Figure 4 and Figure 5 ). This means that the stability and expressibility of the newly designed mRNA sequence can be increased.

[0136] However, by comparing the identity of the CDS regions of the 10 sequences generated by the beam search method and the cluster search method, it can be found that the identity of the CDS regions of the sequences generated by the beam search method is too high (99.8%); on the contrary, without a significant loss of MFE and CAI indicators, the identity of the CDS regions of the sequences generated by the cluster search method (94.6%) is greatly improved (as shown in Tables 1 and Figure 6 As shown, Identity CDS Refers to the CDS regional identity, Average Identity CDS The average degree of identity in the CDS region. If measured by sequence diversity, the sequence diversity generated by the cluster clustering 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 Beam Search and Beam Cluster Search methods (D=100, R=10).

[0138]

[0139] Take the 10 recommended sequences generated by the cluster search method (K=4, D=100, drop rate=0.02) as an example. The average identity of the CDS region of these 10 sequences with the CDS region of the BioNTech vaccine is 87.6%; 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 comparison, the CDS regions of the 10 recommended sequences generated by the beam search method (K=4, D=100) have an average identity of 85.7% with the CDS region of the BioNTech vaccine, and an average identity of 88.6% with the CDS region of the Moderna vaccine; the identity between the sequences is 99.8% (see Figure 8 ).

[0140] In addition, we also studied the effects of the drop rate and search depth D in the algorithm on the diversity of the designed sequences. We found that when the drop rate increased from 1% to 5%, the diversity of the designed sequences decreased. Therefore, choosing a smaller drop rate is beneficial to increase sequence diversity (see Fig. 9 , where Identity CDS Refers to the CDS regional identity, Average Identity CDS It refers to the average identity of the CDS region. However, when the search depth increases from 50 to 500, the diversity of the designed sequences increases slightly, but not significantly (see Fig.10 , where Identity CDS Refers to the CDS regional identity, Average Identity CDS refers to the average identity of the CDS region).

[0141] In summary, the present preferred embodiment relates to a new messenger RNA (mRNA) sequence design method. 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. In this way, its stability can be optimized from the perspective of the overall molecule, and the expressibility of the CDS region can be optimized at the same time. Practice has shown that the mRNA sequence designed using this method has better stability and expressibility than the mRNA sequence designed by the traditional method. This method can be applied to sequence design and optimization in the fields of mRNA therapeutic drugs, mRNA vaccines, gene therapy drugs, etc.

[0142] The various method implementations of the present invention can be implemented in software, hardware, firmware, etc. Regardless of whether the present invention is implemented in software, hardware, or firmware, the instruction code can be stored in any type of computer-accessible memory (e.g., permanent or modifiable, volatile or non-volatile, solid or non-solid, fixed or replaceable media, etc.). Similarly, the memory can be, for example, a programmable array logic (PAL), a random access memory (RAM), a programmable read-only memory (PROM), a read-only memory (ROM), an electrically erasable programmable read-only memory (EEPROM), a magnetic disk, an optical disk, a digital versatile disc (DVD), etc.

[0143] The second embodiment of the present invention relates to a messenger RNA sequence design device. Fig.11 It is a schematic diagram of the structure of the messenger RNA sequence design device.

[0144] Specifically, if Fig.11 As shown, the messenger RNA sequence design device comprises:

[0145] The first sequence generation unit is used to add one codon from all possible codons corresponding to the first amino acid encoding the protein at the end of a given messenger RNA sequence, 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;

[0146] A first judging unit, configured to judge whether the number of sequences in the sequence pool is greater than a search depth D;

[0147] A second sequence generating unit is used for adding one codon from all possible codons corresponding to the next amino acid encoding the protein to the end of each sequence in the sequence pool when the first judging unit judges that the number of sequences in the sequence pool is not greater than the search depth D, thereby generating a group of new sequences; wherein only one codon is added to the end of the sequence, and each codon is added only once;

[0148] The first calculation unit is used to calculate the full-length minimum free energy MFE of each sequence in the sequence pool when the first judgment unit determines that the number of sequences in the sequence pool is greater than the search depth D. whole and codon usage efficiency CAI in the CDS region CDS ;

[0149] A second calculation unit is used to calculate the full-length minimum free energy MFE output by the first calculation unit. whole and the codon usage efficiency CAI of the CDS region CDS , calculating a score for each sequence in the sequence pool, and sorting all sequences in the sequence pool according to the score;

[0150] A filtering unit, used 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 is used 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 is used to judge whether the codons corresponding to all the amino acids of the encoded protein have been added;

[0153] The second sequence generating unit is further configured to, when the second judging unit judges that the codons corresponding to all the amino acids of the encoded protein have not been added, add one codon from all the possible codons corresponding to the next amino acid of the encoded protein to the end of all the sequences in the sequence pool, thereby generating a set of new sequences; wherein only one codon is added to the end of the sequence, and each codon is added only once;

[0154] A 3'UTR sequence adding unit, used for adding 3'UTR sequences to the ends of all sequences in the sequence pool respectively when the second judgment unit judges that the codons corresponding to all the amino acids of the encoded protein have been added;

[0155] The third calculation unit is used to calculate the full-length minimum free energy MFE of each sequence after adding the 3'UTR sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ;

[0156] The fourth calculation unit is used to calculate the full-length minimum free energy MFE output by the third calculation unit. whole and the codon usage efficiency CAI of the CDS region CDS , calculating the score of each sequence after the 3'UTR sequence is added to the sequence pool, and sorting all the sequences after the 3'UTR sequence is added to the sequence pool according to the score;

[0157] The second clustering unit is used to recommend the number R of sequences, cluster all the sequences after the 3'UTR sequence is added to the sequence pool into R categories, and recommend the sequence with the highest score in each category as the designed messenger RNA sequence.

[0158] The embodiment of the present application proposes a messenger RNA (mRNA) sequence design method in the presence of 5'UTR and 3'UTR, which can simultaneously optimize the full-length minimum free energy of an mRNA sequence and the codon usage efficiency of the CDS region, thereby 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 conjunction with the first embodiment. The relevant technical details mentioned in the first embodiment are still valid in this embodiment, and in order to reduce repetition, they are not repeated here. Accordingly, the relevant technical details mentioned in this embodiment can also be applied in the first embodiment.

[0160] It should be noted that the units mentioned in the various device implementations of the present invention are all logical units. Physically, a logical unit can be a physical unit, or a part of a physical unit, or can be implemented as a combination of multiple physical units. The physical implementation of these logical units themselves is not the most important. The combination of 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-mentioned device implementations of the present invention do not introduce units that are not closely related to solving the technical problems proposed by the present invention, which does not mean that there are no other units in the above-mentioned device implementations.

[0161] It should be noted that those skilled in the art should understand that the implementation functions of each unit shown in the implementation methods of the above-mentioned devices can be understood with reference to the relevant description of the above-mentioned corresponding methods. The functions of each unit shown in the implementation methods of the above-mentioned devices can be implemented by a program (executable instruction) running on a processor, or by a specific logic circuit. If the above-mentioned 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 solution of the embodiments of this specification is essentially or the part that contributes to the prior art can be embodied in the form of a software product, which is stored in a storage medium, including several instructions for 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 each embodiment of this specification. The aforementioned storage medium includes: various media that can store program codes, such as a U disk, a mobile hard disk, a read-only memory (ROM, Read Only Memory), a disk or an optical disk. In this way, the embodiments of this specification are not limited to any specific combination of hardware and software.

[0162] Accordingly, the embodiments of the present application also provide a messenger RNA sequence design computing device, which includes a memory for storing computer executable instructions, and a processor; the processor is used to implement the steps in the above-mentioned method implementation methods when executing the computer executable instructions in the memory. Among them, the processor can be a central processing unit (Central Processing Unit, referred to as "CPU"), or other general-purpose processors, digital signal processors (Digital Signal Processor, referred to as "DSP"), Application Specific Integrated Circuit (Application Specific Integrated Circuit, referred to as "ASIC"), etc. The aforementioned memory can be a read-only memory (read-only memory, referred to as "ROM"), a random access memory (random access memory, referred to as "RAM"), a flash memory (Flash), a hard disk or a solid-state hard disk, etc. The steps of the method disclosed in each embodiment of the present invention can be directly embodied as a hardware processor to be executed, or a combination of hardware and software modules in the processor can be executed.

[0163] In addition, an embodiment of the present application further provides a computer storage medium, in which computer executable instructions are stored. When the computer executable instructions are executed by a processor, the various method embodiments of the present application are implemented.

[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 terms "include", "comprise" or any other variants thereof are intended to cover non-exclusive inclusion, so that the process, method, article or device including a series of elements includes not only those elements, but also other elements not explicitly listed, or also includes elements inherent to such process, method, article or device. In the absence of further restrictions, the elements defined by the sentence "including one" do not exclude the existence of other identical elements in the process, method, article or device including the elements. In the application documents of this patent, if it is mentioned that an action is performed according to an element, it means that the action is performed at least according to the element, which includes two situations: performing the action only according to the element, and performing the action according to the element and other elements. Expressions such as multiple, multiple, and multiple include 2, 2 times, 2 kinds, and more than 2, more than 2 times, and more than 2 kinds.

[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 the basis for modification when 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 claimed in this application.

Claims

1. A method for designing a messenger RNA sequence, characterized in that: The following steps are involved: Step S1: adding one codon from all possible codons corresponding to the first amino acid encoding the protein to the end of a given messenger RNA sequence, generating a set of initial sequences, and placing 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; Step S2: determining whether the number of sequences in the sequence pool is greater than the search depth D; If not, execute step S3: add one codon from all possible codons corresponding to the next amino acid encoding the protein to the end of all sequences in the sequence pool, thereby generating 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 added only once; If yes, then execute step S4: calculate the full-length minimum free energy MFE of each sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ; Step S5: According to the full-length minimum free energy MFE whole and the codon usage efficiency CAI of the CDS region CDS , calculating a score for each sequence in the sequence pool, and sorting all sequences in the sequence pool according to the score; Step S6: According to a preset filtering rule or elimination rate, the sequences in the sequence pool that do not meet the requirements are filtered out or the sequences with the lowest ranking in the sequence pool are eliminated, and the remaining sequences are retained; Step S7: clustering the retained sequences into D categories according to the search depth D, retaining the sequence with the highest score in each category in the sequence pool, and eliminating the remaining sequences; Step S8: determining whether all the codons corresponding to the amino acids of the encoded protein have been added; If not, return to step S3; if yes, execute step S9: add 3'UTR sequences to the end of all sequences in the sequence pool; Step S10: Calculate the full-length minimum free energy MFE of each sequence after adding the 3'UTR sequence to the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ; Step S11: According to the full-length minimum free energy MFE whole and the codon usage efficiency CAI of the CDS region CDS , calculating the score of each sequence after the 3'UTR sequence is added to the sequence pool, and sorting all the sequences after the 3'UTR sequence is added to the sequence pool according to the score; Step S12: according to the number R of recommended sequences, all sequences in the sequence pool after adding the 3'UTR sequence are clustered into R categories, and the sequence with the highest score in each category is recommended as the designed messenger RNA sequence.

2. The method according to claim 1, characterized in that The score of the sequence is the minimum free energy MFE of the entire sequence. whole Codon usage efficiency CAI of the CDS region CDS The formula for the score of the 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, characterized in that 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, characterized in that In the step S1, the given messenger RNA sequence is a 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 RNA sequence design device, characterized in that: include: The first sequence generation unit is used to add one codon from all possible codons corresponding to the first amino acid encoding the protein at the end of a given messenger RNA sequence, 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; A first judging unit, configured to judge whether the number of sequences in the sequence pool is greater than a search depth D; A second sequence generating unit is used for adding one codon from all possible codons corresponding to the next amino acid encoding the protein to the end of each sequence in the sequence pool when the first judging unit judges that the number of sequences in the sequence pool is not greater than the search depth D, thereby generating a group of new sequences; wherein only one codon is added to the end of the sequence, and each codon is added only once; The first calculation unit is used to calculate the full-length minimum free energy MFE of each sequence in the sequence pool when the first judgment unit determines that the number of sequences in the sequence pool is greater than the search depth D. whole and codon usage efficiency CAI in the CDS region CDS ; A second calculation unit is used to calculate the full-length minimum free energy MFE output by the first calculation unit. whole and the codon usage efficiency CAI of the CDS region CDS , calculating a score for each sequence in the sequence pool, and sorting all sequences in the sequence pool according to the score; A filtering unit, used 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; A first clustering unit is used 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 judgment unit is used to judge whether the codons corresponding to all the amino acids of the encoded protein have been added; The second sequence generating unit is further configured to, when the second judging unit judges that the codons corresponding to all the amino acids of the encoded protein have not been added, add one codon from all the possible codons corresponding to the next amino acid of the encoded protein to the end of all the sequences in the sequence pool, thereby generating 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 adding unit, used for adding 3'UTR sequences to the ends of all sequences in the sequence pool respectively when the second judgment unit judges that the codons corresponding to all the amino acids of the encoded protein have been added; The third calculation unit is used to calculate the full-length minimum free energy MFE of each sequence after adding the 3'UTR sequence in the sequence pool whole and codon usage efficiency CAI in the CDS region CDS ; The fourth calculation unit is used to calculate the full-length minimum free energy MFE output by the third calculation unit. whole and the codon usage efficiency CAI of the CDS region CDS , calculating the score of each sequence after the 3'UTR sequence is added to the sequence pool, and sorting all the sequences after the 3'UTR sequence is added to the sequence pool according to the score; The second clustering unit is used to recommend the number R of sequences, cluster all the sequences after the 3'UTR sequence is added to the sequence pool into R categories, and recommend the sequence with the highest score in each category as the designed messenger RNA sequence.

9. A messenger RNA sequence design computing device, characterized in that: include: A memory for storing computer executable instructions; as well as, 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 a processor, the steps in the method according to any one of claims 1 to 7 are implemented.

Images