A piRNA Portrait for Blood Cells of Healthy Donors and COVID-19 Patients

1. Introduction

Non-coding RNAs (ncRNA) are functional RNA molecules that are not translated into a protein. There are various types of ncRNAs, including transfer RNAs (tRNAs), microRNAs, small interfering RNAs (siRNAs), PIWI-interacting RNAs (piRNAs), small nucleolar RNAs (snoRNAs), small nuclear RNAs (snRNAs), and long ncRNAs such as Xist and HOTAIR [1]. The functional roles of ncRNAs in gene expression regulation have been widely recognized. They are known to modulate chromatin structure, regulate the assembly and function of nuclear bodies, alter the stability and translation of mRNAs, and interfere with signaling pathways in various ways [2].

The study of ncRNAs is an active area of research, and there is ongoing debate about the functional significance of various ncRNAs [3].

piRNAs are a large subclass of small non-coding RNA molecules [4]. piRNA molecules primarily function in the epigenetic and post-transcriptional silencing of genetic elements in germ line cells by forming complexes with piwi-subfamily Argonaute proteins [5,6]. They are created from loci that function as transposon traps, providing a kind of RNA-based immune system against transposons, which are DNA sequences that can change their position within the genome [6,7,8].

In addition to their role in transposon silencing, piRNAs are also involved in the regulation of other genetic elements in germ line cells. They are primarily expressed in the testes and ovaries of mammals and are essential for germline development and fertility [9]. The dysregulation of piRNAs has been associated with various diseases [10,11], including cancer [12], where they can influence gene expression in somatic cells [13]. It is also known that many piRNAs are dysregulated in tumor tissues, where they can either promote or suppress tumor growth [14]. For instance, piR-54265/PIWIL2 recruits STAT3 and p-SRC through the PIWIL2 PAZ domain to form the PIWIL2/STAT3/p-SRC complex, which facilitates the activation of signal pathways and phosphorylation of STAT3 by p-SRC, thereby promoting tumorigenesis in colorectal adenocarcinoma [15].

To date, there is insufficient data describing the behavior of piRNA in separate cell populations. Despite being previously associated with germinal cell lines [16,17], recent studies show that human piRNA may be related to a wider range of mechanisms. The evidence of an ancestral piRNA machinery that developed into a more complex and sophisticated form in multicellular eukaryotes could be reflected in the existence of PIWI-like machinery in Trypanosomatids and Apicomplexa. According to Horjales et al., ancient PIWI proteins were widely expressed and served purposes apart from germline preservation [18].

In a recent study of peripheral leukocytes of rheumatoid arthritis patients it has been shown that two immunoregulation piRNAs (piR-hsa-27620 and piR-hsa-27124) were significantly elevated in patients with rheumatoid arthritis [19]. Latest advancements in sequencing and molecular technologies, along with in silico studies have improved the annotation and functional validation of ncRNAs, enabling a better understanding of their roles in specific cell types and diseases [20]. However, the functional significance of ncRNAs, especially in poorly conserved species, remains a subject of ongoing research.

Given the fact how significantly COVID-19 alters the biochemical and morphologic features of blood cells [21,22,23], and taking in account the promising research on interrelations between ncRNA immune regulatory features and COVID-19 [24,25,26], it is necessary to state the importance of studying piRNA expression in blood cells of SARS-CoV-2 infected patients. The aim of this research was to analyze and compare the expression of piRNA in sorted blood cells in healthy control donors and in COVID-19 patients.

2. Material and methods

2.1. Patients and Data Collection

Eleven research participants in total were split up into two groups: six healthy donors and five COVID-19 patients. The six healthy donors in our control group had a mean age of 51 ± 22. The five COVID-19 patients, with an average age of 67 ± 11, arrived at our institution between April 1, 2022, and August 23, 2022. Real-time RT-PCR was used in the laboratory to confirm that all patients had SARS-CoV-2 infection. Two patients were transferred to the infectious disease unit and three patients were admitted to the intensive care unit. The following criteria were met by severe patients for ICU admission: body temperature ≥ 39 °C, respiration rate ≥ 30/min, and oxygen saturation (SpO2) ≤ 93%. These criteria are outlined in our national clinical recommendations.

2.2. Cell sorting and scanning electron microscopy

Sample staining for flow cytometry and subsequent sorting was carried out according to a standardized technique [23].

Erythrocyte sorting was carried out by the following technique:

10 µl whole blood was added to 100 µl PBS, 2 µl anti cd235 (IM2212U) was added, 2 µl cd41 (A07781), cd45 (A79392). 20 minutes incubated in dark, dissolved to 1 ml phosphate-buffered saline (PBS), next, fluorescent sorting was carried out up to 5 million events on a MoFlo Astrios EQ flow cytometer sorter (equipped with 405, 488, 645).

Leukocyte sorting was carried out by 2 different methods as follows:

100 µl of heparinized peripheral blood was taken, 1 ml of VersaLyse lysing solution (Beckman Coulter, USA) was added to each tube, stirred and incubated for 10 min. 10 µl of appropriate monoclonal antibodies were added: CD45-APC-AlexaFluor®750+ (A79392 Beckman Coulter, USA), CD16-PC7+(6607118 Beckman Coulter, USA), CD14-PC5.5+(A70204 Beckman Coulter, USA) was used to distinguish lymphocytes and monocytes.

The samples were stirred by vortex mixing and incubated for 15-20 minutes at room temperature in the dark. Afterwards, 1 mL of phosphate buffer (PBS) was added and centrifuged for 5 min at 1500 rpm, the supernatant was removed and 1 mL of buffer solution was added to each tube. Resuspended precipitant was analyzed on a MoFlo Astrios EQ flow cytometer sorter (equipped with 405, 488, 645).

For isolating subpopulations of granulocytes - eosinophils, basophils and neutrophils the duraclone IM granulocyte antibody panel (B88651) (Beckman Coulter, USA) was used.

Sorting was performed through 70 μm diameter nozzles, 50 000-100 000 events for the leukocyte cell population, into sterile 12x75 eppendorfs.

A standardized protocol for sample preparation for scanning electron microscopy (SEM) was followed [27]. Further analysis was performed using a Zeiss Merlin scanning electron microscope in high resolution mode with EHT 0.400 kV and a magnification of 1.00 KX.

2.3. RNA separation, library preparation and next generation sequencing

RNA was isolated\extracted using ExtractRNA reagent (Evrogen) according to the provided protocol. After extraction, RNA was dissolved in 10 µl RNAse free water. The quality of leukocyte RNA was checked on a TapeStation instrument (Agillent). Only samples with RIN>5 were used for the preparation of leukocyte libraries. 15-90 ng of RNA was taken in each reaction. Reversion, preparation of short RNA libraries, and circularization were performed using the Small RNA Library Prep Kit (BGI, 1000006383) according to the manufacturer’s recommendations. Libraries were purified by electrophoresis in a polyacrylamide 6% gel, cutting out the target band corresponding to the length of the target fragment 18-50b (library length 105-140b). Sequencing was performed on a DNBSEQ-G400 (BGI) instrument using the DNBSEQ-G400RS High-throughput Sequencing Set （Small RNA FCL SE50*) (BGI, 1000016998). After sequencing, the amount of information obtained from each sample was 400-1200MB.

2.4. Bioinformatics and statistical analysis

Bioinformatics analysis was performed in miRMaster 2 version 2.0.0 [28], a comprehensive analysis framework for small RNA-seq data. Current version of miRmaster used Ensembl [29], RNACentral piRNA [30], NCBI RefSeq [31], ENA [32], GeneCards [33], PirBase [34] databases for analyzing RNAseq FASTQ data. Statistical analysis was performed in RStudio (version 2023.09.1+494.pro2) [35]. Reads per million (RPM) was used as a normalized value of expression. Statistical analysis for the results was executed by applying the Wilcoxon–Mann–Whitney test. A p-value < 0.05 was considered statistically significant. For RNAseq volcano plot the p-value was log transformed to log10 (1.31) for data presentation. Heatmap raw expression data was preliminarily log transformed. Data visualization, images, and charts were made with RStudio ggplot2 [36] and tidyverse [37] open-source collection of packages. All raw data used for statistical analysis and visualization is presented in the supplementary materials.

3. Results

3.1. piRNA expression of sorted healthy donor blood cells

In order to compare piRNA expression in blood cells we sorted them by distinguishing the following cell types: erythrocytes (n=5), monocytes (n=5), lymphocytes (n=4), eosinophils (n=2), basophils (n=2), neutrophils (n=6). The cell sorting procedure is described in materials and methods. For accuracy, we checked the quality of sorting with SEM. As can be seen in the presented images (Figure 1 panel A) morphological features of cells correspond with sorted cell types [22,23]. We then isolated RNA from sorted cells, prepared small RNA libraries, performed NGS of small RNA and subjected the raw RNAseq data to bioinformatic analysis in miRMaster 2, details are described in materials and methods. All annotated piRNA data is presented in the supplementary materials. We compared the major RNAs overexpressed in different blood cells of healthy control donors. The results are presented in Figure 1 panel A (piechart). Additionally, log transformed raw expression data is presented in the heatmap (Figure 1 panel B), where 5 (red) indicates maximum expression values, and <1 (green) minimum. 4 of the most abundant piRNA were expressed relatively equally in all cell types, they included piR-49145, piR-49144, piR-49143 and piR-33151.

3.2. piRNA differential expression in blood cells of patients with COVID-19

Next, we analyzed the piRNA differential expression in healthy donors and patients with severe COVID-19 in three distinguished sorted cell types erythrocytes (n=5), monocytes (n=5) and lymphocytes (n=4), the results of the comparison are presented as a volcano plot (Figure 2) where positive increasing log2 fold change corresponds with gene downregulation, and negative values correspond with gene upregulation. All raw data used for analysis is presented in the supplementary materials section.

By reviewing piRNA expression in each group of sorted cells from severe COVID-19 patients we observed 3 upregulated and 10 downregulated piRNAs in erythrocytes (Table 1). Monocytes were presented with a larger amount of statistically significant piRNA, 4 upregulated and 35 downregulated (Table 2). In lymphocytes, each piRNA was upregulated (Table 3).