{'fraction_non_alphanumeric': 0.03733749874030031, 'fraction_numerical': 0.03139171621485438, 'mean_word_length': 4.70316091954023, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 1, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Mutations in the nonstructural proteins of SARS-CoV-2 may contribute to adverse clinical outcome in patients with COVID-19

Background From late March through April 2021, we experienced a cluster of patients with COVID-19, named “Cluster K”, with rapid severe illness compared with those who were infected before. Methods Patients with COVID-19 who were enrolled in this study were divided into two groups: 66 patients from November 2020 to March 2021 (group A) and 37 patients whose infection links were traced from Cluster K (group B). The primary outcome was mortality rate, and the secondary outcome was maximal oxygen flow rate as the severity of the disease. Viral genome sequences were compared between the two groups. Results Mortality rates were 6.1% in group A and 16.2% in group B (odds ratio: 2.97, 95% confidence interval: 0.65–15.38). The patients in group B required high oxygen flow rate (O2 ≥10 l/min) in the earlier clinical course (P = 0.029). Viral genome sequences revealed five amino acid mutations; of these, four were found on three nonstructural proteins (NSPs): one in nsp3 and nsp15, two in nsp6 (one of them is near the potential sites under positive selective pressure). Another one was on the S protein. Conclusion This study suggests that mutations in NSPs, especially nsp6, are associated with adverse clinical outcome in patients with COVID-19.


Introduction
COVID-19, which is caused by SARS-CoV-2, initially broke out in Wuhan, China on December 2019 and has subsequently spread around the world ( Gorbalenya et al., 2020 ). SARS-CoV-2 is an enveloped virus with a positive-sense single-stranded RNA genome of approximately 30 kb. Two-thirds of the viral genome at the 5 terminus contains ORF1a and ORF1b, which encodes 16 nonstructural proteins (NSPs). These NSPs play crucial roles in viral replication and evasion of host immune systems ( Maier et al., 2015 ). The other one-third of the 3 end contains structural genes such as spike glycoprotein (S protein), envelope protein, mem-brane protein, nucleocapsid protein, and several accessory proteins ( Hagemeijer et al., 2010 ). Because replication of RNA viruses typically has a low fidelity, the viral genome of SARS-CoV-2 has accumulated mutations at an average of two nucleotides per month ( Duchene et al., 2020 ). On February 25, 2021, the World Health Organization (WHO) defined a variant of interest (VOI) as an isolate that is detected in several countries that changes phenotype under certain conditions and defined a variant of concern (VOC) as an isolate in VOIs if it has been demonstrated to be associated with an increase in transmissibility or virulence ( World Health Organization, 2021 ). For example, VOC-202012/01 (Pango lineage B.1.1.7), called the "alpha variant," has a higher transmissibility and mortality than a reference isolate ( Challen et al., 2021 ;Davies et al., 2021 ). These variants, however, were not focused on mutations in NSPs but in structural proteins, especially the S protein. As of May 2022, COVID-19 is responsible for 8.8 million cases and approximately 30,0 0 0 deaths in Japan ( Ministry of Health, Labour and Welfare, 2022a ) ( https://www.mhlw.go.jp/stf/ covid-19/kokunainohasseijoukyou.html ). The most remarkable clinical feature of this disease is its heterogeneity in clinical manifestations, ranging from no symptoms to critical illness ( Grasselli et al., 2020 ). Several risk factors for symptom severity have been known, such as age (over 65 years old), male sex, obesity, history of smoking, and comorbidities (including hypertension, diabetes mellitus, respiratory disease, cardiovascular disease, chronic kidney disease, and cancer) ( Gansevoort and Hilbrands, 2020 ;Liang et al., 2020 ;Popkin et al., 2020 ;Zheng et al., 2020 ). In Japan, several treatment options are currently recommended on the basis of symptom severity: molnupiravir, nirmatrelvir and ritonavir, casirivimab and imdevimab, sotrovimab, remdesivir, dexamethasone, and baricitinib ( Ministry of Health, Labour and Welfare, 2022b ). In clinical practice, methylprednisolone is also used for patients who are critically ill when they have a poor response to dexamethasone therapy or when baricitinib is not available ( Pinzón et al., 2021 ).
In Asahikawa City, Hokkaido Prefecture, located in northern part of Japan, a cluster of COVID-19 cases occurred in several bars where a group of customers enjoyed "karaoke" during daytime from the end of late March through April 2021 (referred to below as Cluster K), and 66 cases were involved into the infection according to the city's survey (Asahikawa Health Center, 2021) . We found that the progression of disease among the patients in this cluster was faster compared with patients who were infected before that. Therefore, to identify factors associated with the progression of disease, we examined the viral genome sequences along with the clinical characteristics and treatment strategies, and then compared those with other cases.

Participants
All adult patients with COVID-19 who were admitted to Asahikawa City Hospital from November 2020 to April 2021 were enrolled in this retrospective study. All positive results were confirmed by real-time reverse transcription polymerase chain reaction (RT-PCR) for the presence of SARS-CoV-2 in saliva or nasopharyngeal swab samples. The recovered patients were discharged from the hospital following the discharge criteria provided by Japanese Ministry of Health, Labour and Welfare when two criteria have been met: 10 days passed from the symptom onset and 72 hours passed from the symptom resolution or when two consecutive negative PCR results are confirmed.

Ethics
This study was approved by the ethics committee of Asahikawa City Hospital (approval number 5, 2021). The opt-out method was used to obtain patient consent in this study. We provided the patients with information explaining the proposed research plan through the website of Asahikawa City Hospital.

Clinical data collection and definitions
The epidemiological data, medical history, underlying comorbidities, symptoms and signs at admission, oxygen flow rate, treatment, and clinical outcomes were obtained from electronic medical records. Because the date of onset was unclear for some patients who were asymptomatic, we considered the date on which the positive PCR result was obtained as day 0. Clinical outcomes were followed up until discharge or death. The primary outcome was mortality rate, and the secondary outcome was severity of the disease, stratified by maximal oxygen flow rate into O 2 ≥1 l/min, ≥5 l/min, and ≥10 l/min. The days to reach its oxygen dose was calculated as occurrence of an event, and cumulative probability was plotted using the Kaplan-Meier method.

Viral genome sequence
The stored RNA extracts were used as clinical samples. Fiststrand cDNA was synthesized by using a PrimeScript IV first strand cDNA Synthesis Mix (TaKaRa Bio) with random hexamer primers and extracted RNA, according to the manufacture's protocols. The whole-genome sequencing of SARS-CoV-2 was carried out as reported previously ( Torii et al., 2021 ). Briefly, a total of nine gene fragments shown in Figure S1 were amplified with synthesized cDNA, specific primer sets for SARS-CoV-2 and PrimeSTAR GXL DNA polymerase (TaKaRa Bio). Then, the amplified products were directly sequenced in both directions by using the ABI PRISM 3130 Genetic Analyzer (Applied Biosystems) with specific primers. The primer sets used in this study are listed in Table S1. Our raw sequence data were converted to FASTQ format by in-house Python script with BioPython module ( Cock et al., 2009 ), followed by trimming first 50 bps and the region after 900 bps by FASTX-Toolkit ( http://hannonlab.cshl.edu/fastx _ toolkit/ ). Trimmed sequence reads were mapped to the SARS-CoV-2 reference genome (accession number MN908947.3).
All SARS-CoV-2 genome sequence data, including our assemblies, were mapped to the reference genome (accession number MN908947.3). On the basis of the SAM format mapping output, SNPs and short INDELs were detected and categorized into synonymous, nonsynonymous, or intergenic mutations by in-house Perl scripts, which simply picked up the differences from the reference sequence for each genome. For the phylogenetic analysis, a typical genome sequence of each variant, named alpha, beta, gamma, delta, epsilon, zeta, eta, theta, iota, kappa, and lambda by WHO, were chosen from the published data on the basis of the following criteria: neither 'N' nor other mixed bases were included in open reading frames (ORFs), number of substitution sites in ORFs were less than 100, all nonsynonymous mutations that characterized each variant were included, and all other substitutes observed in more than 95% of the samples that satisfied the above conditions were also included. The nonsynonymous mutations that characterized each variant are shown in Table S2, created on the basis of the mutation prevalence information from outbreak.info ( outbreak.info, 2020 ). ORF10 was excluded from the analyses because one of our assemblies did not completely cover that region. On the basis of the alignment, suitable nucleotide substitution models for each partition were selected by ModelTest-NG ( Darriba et al., 2020 ), followed by maximum likelihood phylogenetic analyses using the selected substitution models by RAxML-NG ( Kozlov et al., 2019 ). The constructed phylogenetic tree was drawn and edited by MEGA X ( Kumar et al., 2018 ).

Statistical analysis
We performed group comparisons using the Fisher's exact test for categorical variables and independent group t -test for continuous variables. The categorical variables were expressed as frequencies (percentage) and the continuous variables as the median with interquartile range. Odds ratios (ORs) with 95% confidence interval (CI) were calculated for the rates of mortality and the required oxygen flow rate. Cumulative probabilities of oxygen therapy were described using the Kaplan-Meier analysis and compared using the log-rank test. P -values of 0.05 or less were considered statistically significant. All statistical analyses were performed with EZR (Saitama Medical Center, Jichi Medical University, Saitama, Japan, version 1.54) ( Kanda, 2013 ), which is a graphical user interface for the R software (The R Foundation for Statistical Computing, Vienna, Austria).

Results
A total of 191 patients diagnosed with COVID-19 were admitted to our hospital from November 2020 to April 2021. The patients were divided into two groups: group A, the patients who were infected from November 2020 to March 2021, and group B, the patients whose infection links were traced from Cluster K. We included the patients whose age ranged between 65 and 89 years for the following reasons: (i) patients who are younger than 65 years old are not necessarily requested for hospitalization in our city's health care system and (ii) in group A, many patients over 90 years old were transferred from a convalescent hospital where COVID-19 outbreak had occurred; thus, the age profiles would be much different from group B if the patients aged over 90 years old were included in this group. Finally, 66 patients were assigned to group A and 37 patients to group B ( Figure 1 ).
The clinical characteristics and treatment strategies used during hospitalization are shown in Table 1 . Although there were no significant differences in clinical characteristics between the two groups, obesity and lifestyle-related diseases such as hypertension, dyslipidemia, and diabetes tended to be more prevalent in group B. In treatment strategies, remdesivir tended to be used more frequently in group A and methylprednisolone in group B.
A total of four of the 66 patients (6.1%) in group A and six of the 37 patients (16.2%) in group B died of COVID-19 (OR: 2.97, 95% CI: 0.65-15.38) ( Table 2 ). When limited to males, the differences in number of deaths became larger (6.3% vs 29.4%, OR: 5.99, 95% CI: 0.84-71.2). The clinical characteristics and the treatment of the ten deceased patients are shown in Table S3. Because no significant difference was found in mortality rate, the severity of disease,   stratified by maximal oxygen flow rate was evaluated as a secondary outcome. The cumulative probability of the patients who required low dose oxygen therapy (O2 ≥1 l/min) was similar in both groups ( Table 3 , Figure 2 A). However, the patients in group B were likely to require more oxygen therapy in the earlier course of hospitalization ( Figure 2 B, 2 C). For high-dose oxygen therapy (O2 ≥10 l/min), a significant difference was observed ( P = 0.029). With mortality, males in group B showed worse outcomes in the severity of the disease ( Table 3 ).
Assuming that these differences in clinical outcome would be due to viral factors, we examined the whole viral genome sequences, representing two samples in each group, named Asahikawa_0108 and Asahikawa_0122 isolated in January 2021 and Asahikawa_0404 and Asahikawa_0417 isolated in April 2021. The mutation sites, which differed in these four samples, are summarized in Table 4 . Focusing on the nonsynonymous mutations, five amino acid changes were identified (represented in bold). Of these, four amino acid changes were found on NSPs and one was on the S protein, as illustrated in the viral genome scheme of SARS-CoV-2 ( Figure 3 ). To assess strain relatedness, the phylogenetic analysis was carried out. The isolates in this study were closely related to the strain previously isolated in Japan (accession number  G  A  A  A  D1323G  6380  C  T  C  C  L2039F  6433  T  C  C  C  Syn  6532  G  G  T  T  E2089D  6997  T  C  C  C  Syn  8371  T  G  G  G  Q2702H  8917  C  T  T  T  Syn  9207 T Orf1b Position number based on alignment to NCBI reference sequence: NC_045512 (Wuhan-Hu-1). Amino acid changes identified between two groups were represented in bold.

Figure 2.
Cumulative probability of the patients who required oxygen therapy. The days to reach its oxygen dose was calculated as occurrence of an event, and cumulative probability was plotted using the Kaplan-Meier method. The log-rank test was used for statistical analysis.  BS0 0 0756.1; Pango lineage B.1.1.214) but not to any of the variants labeled by WHO ( Figure 4 ).

Discussion
In this study, we experienced the cases of a single COVID-19 cluster, in which the conventional treatment was unsuccessful, and the hospitalized patients became more severely ill; thus, we performed the retrospective analysis to find out the causative factors. The patients in group B had a higher mortality rate than those in group A, but no significant difference was indicated. For the severity of disease, the patients in group B required a higher oxygen flow rate in early course of admission. Furthermore, the differences became larger when limited to males, with 41.2% of male severe patients in group B against 6.3% in group A.
In terms of treatment strategies, we were likely to use methylprednisolone monotherapy for severe cases in group B. This was because we have found, in the winter of 2020, that patients who received methylprednisolone therapy appeared to have a better response than those who were treated with remdesivir after poor response for dexamethasone. In this study, similar results were obtained when the secondary outcome was compared only for patients without remdesivir (data not shown). Therefore, the viral factors rather than the difference of the treatments could have an impact on the clinical outcome.
Whole viral genome sequences revealed five nonsynonymous mutations by comparison of the isolates with each group. Interestingly, four of them were found on NSPs. The sole structural mutation, D1153G on the S protein, was located near the C-terminus and outside the range of three-dimensional structures obtained by crystallography, suggesting neither the head nor stalk of the S protein. Therefore, D1153G is not expected to be involved in the virulence of the disease. Of the four mutation sites in NSPs, two were found on nsp6 (T3598I and Y3722C), and the rest were on nsp3 (E2089D) and nsp15 (A2143S). To our knowledge, none of these mutations are known to be the factors that determine severity of the disease in clinical practice, and it is unclear whether the mutations affect the functions of these NSPs. However, all of them have been reported to play roles in the suppression of type 1 interferon in host cells ( Shemesh et al., 2021 ;Shin et al., 2020 ;Xia et al., 2020 ;Yuen et al., 2020 ). In addition, nps6 is involved in the formation of autophagosome in host cells. The potential sites under positive selective pressure have been found on nsp6 near T3598I according to an evolutionary analysis on SARS-CoV-2 genome sequences of 351 clinical samples ( Benvenuto et al., 2020 ). Therefore, the mutations in nsp3, nsp6, and nsp15 may be responsible for its function to interact with host immunity and autophagy. In partic-ular, nsp6 is known to interact with the sigma-1 receptor, which is a transmembrane endoplasmic reticulum protein that modulates activity of multiple effector proteins and is a potential target for antiviral drugs ( Gordon et al., 2020a ). Hydroxychloroquine has been reported to inhibit viral activity through the interaction with the sigma-1 receptor ( Gordon et al., 2020b ), suggesting that nsp6 could play an important role clinically. Further investigations are required to reveal their functions and involvement in virulence.
The phylogenetic analysis revealed that the isolates in this study are derived from Pango lineage B.1.1.214, which was the dominant strain in the so-called "third wave," the epidemic period from October to December 2020 in Japan ( National Institute of Infectious Diseases, 2021 ). Considering the fact that any strains identical to the viral sequence of group B cannot be found in more than 1,0 0 0,0 0 0 registered sequences in NCBI, it is likely that the viruses evolved locally from January to April 2021 in Hokkaido. Since May 2021, alpha variant has been dominant throughout Japan, including Hokkaido; 78% of the samples isolated from April 26 to May 2, 2021, in Hokkaido were already the alpha variant ( Ministry of Health, Labour and Welfare, 2021 ). Therefore, the B.1.1.214 strains in Asahikawa seem to have been eliminated. This is similar to the situation in influenza A virus, for example, where pandemic influenza A (H1N1) pdm09 replaced Russian flu (H1N1) since 2009 ( Blyth et al., 2010 ;Broor et al., 2012 ).
There are some limitations to this retrospective study. First, the recommended therapy for COVID-19 has been changing dynamically. Because the best treatment strategies are selected on the basis of the evidence at that time, it is difficult to compare clinical outcome between different periods of time in real-life practice. Second, this is a single-institution study; not all cases in Cluster K have been assessed. The mortality and severity rate would be different when including the cases admitted to other hospitals. Third, we examined only two samples from each group in this study. Due to the limited human and financial resources, we were unbale to perform the viral genome analysis of more samples.
In conclusion, this study demonstrated that some viral factors, four nonsynonymous mutations in nsp3, nsp6, and nsp15 but not on the S protein, may contribute to adverse clinical outcome in patients with COVID-19. Although scientists and physicians all over the world have focused on the S-protein mutations, to comprehensively understand the epidemiology of COVID-19, we should have an insight into other types of mutations such as NSPs.

Data availability
All data used and/or analyzed during this study are available from the corresponding author on reasonable request.