Frequently Asked Questions on RNA-seq

Inquiry

In recent years, the significant reduction in sequencing costs and the burgeoning interest in transcriptomics research have established RNA-seq as the preferred technique for advancing projects in various molecular biology research groups. During our participation in transcriptome research projects, we have observed varying levels of uncertainty and confusion among researchers regarding the analysis results. To address these concerns, we have compiled a series of common questions related to transcriptome projects for researchers' reference. Researchers can select questions of interest or those needing clarification based on their individual needs. Should any additional questions arise during the reading, we encourage reaching out for free consultation, where we will provide expert assistance and guidance.

1. How is the consistency of biological replicates determined? Methods and formulas for biological replicate statistics

Pearson correlation coefficient (r) is a key metric for evaluating the consistency of biological replicates. Ideally, the correlation coefficient (r²) should be ≥0.92 for high-quality biological replicate experiments. Given the influence of factors such as individual variability, sampling environment, timing, and operator expertise on sequencing data, a correlation coefficient (r²) ≥0.8 is generally considered acceptable.

Pearson correlation coefficient: Also known as the product-moment correlation coefficient, it was introduced by the British statistician Karl Pearson in the 20th century as a measure of the linear correlation between two variables.

2. Should DEGs be represented by Transcripts or Unigenes?

3. What does transcript-id signify? Why might a gene have multiple transcript-ids?

4. How are reads potentially identified as miRNA calculated in miRNA identification? What factors influence the mrd value? How are miRNA isomers in different tissues handled?

miRNA identification involves aligning sequences with databases such as Rfam, miRBase, RepBase, and Exon/Intron databases to obtain small RNA annotation information, serving as the foundation for predicting novel miRNAs. The miRDeep2 software is utilized for identifying both known and novel (conserved and non-conserved) miRNAs. MiRDeep2 extends the positions of reads mapped to the genome by 75 and 15 bp at both ends, respectively, for structural prediction. The software differentiates between potential and highly probable miRNAs based on the mrd value, where mrd > -10 indicates potential miRNAs, and mrd > 0 indicates highly probable miRNAs. The distribution of reads on the genome and the free energy of base pairing influence the mrd value.

5. How are differentially expressed genes calculated in projects with biological replicates?

6. How is the proportion of exons, introns, and intergenic regions evaluated to assess library quality?

7. Factors influencing Contig assembly results

Species specificity: Unique characteristics of the species can affect the assembly process.

Sequencing quality: High-quality sequencing data yield better assembly results.

Data volume: Greater sequencing data volume increases the number of short Contigs due to a higher SNP heterozygosity rate. The Trinity assembly strategy constructs a K-mer library from Reads, selects the most frequent K-mers, and extends them using k-1 overlap. The more reads sequenced, the higher the SNP heterozygosity rate, resulting in more short Contigs after extension.

Assembly parameters: The choice of parameters for filtering low-frequency K-mers affects the number of filtered K-mers and, consequently, the number of remaining short Contigs. Thus, even with the same software (e.g., Trinity), differences in assembly parameters can lead to incomparable results.

Assembly quality: The quality of assembly is primarily assessed by the number and N50 of assembled Unigenes. Generally, if the number of assembled Unigenes falls within a reasonable range (e.g., within 100,000), a higher N50 indicates better assembly results.

8. Differences Between Contig and Transcript in Transcriptome Sequencing

Transcriptome sequencing generates raw data consisting of numerous reads. During sequence assembly, overlapping reads are combined into larger fragments known as contigs. These contigs are further aligned, and paired-end reads can ascertain the association of different contigs derived from the same transcript and their relative distances. This assembly process culminates in sequences that cannot be further extended at both ends, termed unigenes. A transcript refers to the actual RNA molecule expressed from the genome.

9. Do Different ID Numbers Represent the Same Gene? Why Do Different IDs Have the Same Functional Annotation?

Distinct ID numbers generally correspond to different genes. However, genes annotated with identical functions may arise from two main scenarios: firstly, some long genes may not be fully assembled, resulting in multiple small fragments that annotate to the same functional protein; secondly, although the nucleotide sequences of genes may differ, their protein sequences can exhibit significant similarity or contain similar functional domains, leading to the same functional annotation when aligned with reference protein sequences.

10. Why Do Multiple Unigenes Have the Same Annotation Despite Different Sequence Lengths and Low Similarity?

11. Interpretation of Fields: transcript_id, gene_id, length, effective_length, expected_count, TPM, FPKM, IsoPct

A single unigene may correspond to multiple transcripts. The fields are defined as follows:

transcript_id: The identifier for the assembled transcript.

gene_id: The identifier for the unigene.

length: The length of the unigene.

effective_length: The average length of all transcripts corresponding to a gene.

TPM: Transcripts per million, calculated as the number of reads mapping to a unigene multiplied by 10^6, divided by the total number of reads.

FPKM (Fragments Per Kilobase of transcript per Million mapped reads): Equivalent to RPKM (Reads Per Kilobase of transcript per Million mapped reads) for paired-end data, calculated as the number of fragments mapping to the transcript divided by the total number of fragments mapping to all transcripts, normalized by transcript length.

IsoPct: The proportion of expression of a particular transcript relative to the total expression of all transcripts from the same gene.

12. Selection of Nucleotide Sequence for a Given ID with Multiple Sequences

13. Methodology for Differential Expression Threshold and Differences Between p-value and FDR

In transcriptome analysis, the Benjamini-Hochberg method is employed to correct the raw p-values obtained from hypothesis testing. The corrected p-values, known as the False Discovery Rate (FDR), are used as a key metric for selecting differentially expressed genes. This correction method mitigates the false positive rate associated with conducting multiple independent statistical tests on gene expression levels. The FDR is not directly convertible to the p-value and is calculated using a program written in the R language under the Linux operating system. For screening, an FDR < 0.01 and a fold change (FC) ≥ 2 are set as the selection criteria.

14. Calculation of Fold Change for Genes with Zero Expression

15. Definition of Known miRNA, Conserved miRNA, and Newly Predicted miRNA

16. Identical Precursor and Mature Sequences with Different Names: Causes

This phenomenon occurs due to differences in chromosomal locations. For instance, in the miRBase database, the IDs hsa-mir-1233-1 and hsa-mir-1233-2 represent miRNAs with identical precursor and mature sequences, but located at different genomic loci. As a result, software distinguishes them as separate small RNAs, even though the reads aligned to them and their expression levels are the same across samples.

Example:

Precursor Sequences:

hsa-mir-1233-1: GUGAGUGGGAGGCCAGGGCACGGCAGGGGGAGCUGCAGGGCUAUGGGAGGGGCCCCAGCGUCUGAGCCCUGUCCUCCCGCAG

hsa-mir-1233-2: GUGAGUGGGAGGCCAGGGCACGGCAGGGGGAGCUGCAGGGCUAUGGGAGGGGCCCCAGCGUCUGAGCCCUGUCCUCCCGCAG

Mature Sequences:

hsa-mir-1233-1 and hsa-mir-1233-2, miR-1233-5p: AGUGGGAGGCCAGGGCACGGCA

hsa-mir-1233-1 and hsa-mir-1233-2, miR-1233-3p: UGAGCCCUGUCCUCCCGCAG

Chromosomal positions:

hsa-mir-1233-1: chr15, 34382069-34382150

hsa-mir-1233-2: chr15, 34528290-34528371

17. Distinguishing Known and Novel lncRNAs in Sequencing Data

For well-characterized species, such as humans, rats, and mice, which possess comprehensive known lncRNA databases, the identification of known and novel lncRNAs can be accomplished as follows:

Positional Relationship: Known lncRNAs are identified based on their specific genomic locations. If a predicted lncRNA overlaps with a known lncRNA, it is considered the same.

Sequence Comparison: The predicted lncRNA sequences are aligned with known lncRNA sequences in databases. If the sequence similarity meets a certain threshold, the predicted lncRNA is identified as a known lncRNA. The NONCODE database, which contains lncRNAs from various species, serves as a reference for this analysis.

18. Discrepancies Between Sequencing Results and qPCR Validation

To reconcile discrepancies between RNA-seq data and qPCR validation results, consider the following factors:

Sample Consistency: Verify whether the samples used in validation are from the same batch and if their up/down-regulation trends align with the sequencing results (e.g., T01 vs. T02, where T01 is the control and T02 is the experimental group).

Expression Levels and Contamination: Ensure that the expression levels of the target genes are not excessively low, which could affect detection accuracy. Also, confirm that neither the samples nor the primers used in the qPCR experiments are contaminated, as this could compromise the results.

Accuracy of qPCR Results: If discrepancies persist after addressing the above factors, re-examine the accuracy of the qPCR results. A single inconsistent gene may not indicate a problem; however, if multiple genes show inconsistent results, it may be necessary to reassess the validity of the sequencing data and the qPCR results.

19. Converting SAR Format Data from NCBI to FASTQ Format

20. Small Peaks in Transcriptome Insert Size Detection

Small peaks appearing during the detection of insert sizes in reference-based transcriptome analysis are considered normal. These peaks arise when paired-end reads map to two adjacent exonic regions on the reference genome, resulting in a measured insert size longer than the actual sequence. This phenomenon is expected in reference-based analyses, where small peaks to the right of the main peak are normal. However, such peaks are considered anomalous in other types of analyses or if they appear to the left of the main peak.

* For Research Use Only. Not for use in diagnostic procedures.