Ute RNA-Seq Data Methods

PCA plot untreated/wt vs parental cells; 1x für WaGa cell line und 1x für MKL-1 cells

Heatmap untreated/wt vs parental; 1x für WaGa cell line und 1x für MKL-1 cells

Volcano plot untreated/wt vs parental; 1x für WaGa cell line und 1x für MKL-1 cells

Distribution of different RNA Species untreated/wt and parental; 1x für WaGa cell line und 1x für MKL-1 cells

RNA fragmentation patterns: is EV-RNA full length or fragmented?

RNA binding protein motifs: do we find specific motifs in EV-RNA?

miRNA target analysis

Input files (use R 4.3.3)

    ~/DATA/Data_Ute/Data_RNA-Seq_MKL-1+WaGa/merged_gene_counts_40samples.txt
    ~/DATA/Data_Ute/Data_RNA-Seq_MKL-1+WaGa/merged_gene_counts_virus_rounded.txt

Common processing for the data MKL+1 + WaGa

    #[1] "/home/jhuang/mambaforge/envs/r_env/lib/R/library"

    setwd("/home/jhuang/DATA/Data_Ute/Data_RNA-Seq_MKL-1_WaGa/results_2025_1/")

    library("AnnotationDbi")
    library("clusterProfiler")
    library("ReactomePA")
    #library("org.Mm.eg.db")
    library(DESeq2)
    library(gplots)

    d.raw_human<- read.delim2("../merged_gene_counts_40samples.txt",sep="\t", header=TRUE, row.names=1)
    colnames(d.raw_human)<- c("gene_name","X042_MKL.1_wt_EV","MKL.1_RNA_147","X042_MKL.1_sT_DMSO","MKL.1_RNA","X0505_MKL.1_scr_DMSO_EV","X0505_MKL.1_sT_DMSO_EV","MKL.1_EV.RNA_87","MKL.1_EV.RNA_27","X042_MKL.1_scr_Dox_EV","X042_MKL.1_scr_DMSO_EV","MKL.1_EV.RNA_118","X0505_MKL.1_scr_Dox_EV","MKL.1_EV.RNA","X042_MKL.1_sT_Dox","X0505_MKL.1_sT_Dox_EV","MKL.1_RNA_118","MKL.1_EV.RNA_2",      "Geneid.1","gene_name.1","X1605_WaGa_sT_DMSO_EV","WaGa_EV.RNA_118","WaGa_EV.RNA","X2706_WaGa_scr_Dox_EV","X2706_WaGa_sT_DMSO_EV","X1107_WaGa_wt_EV","X1107_WaGa_sT_Dox_EV","WaGa_RNA","X1605_WaGa_scr_DMSO_EV","X2706_WaGa_scr_DMSO_EV","WaGa_EV.RNA_226","X2706_WaGa_sT_Dox_EV","X1605_WaGa_scr_Dox_EV","X1605_WaGa_wt_EV","X1605_WaGa_sT_Dox_EV","X1107_WaGa_scr_DMSO_EV","WaGa_EV.RNA_2","WaGa_EV.RNA_147","WaGa_RNA_118","X1107_WaGa_scr_Dox_EV","X1107_WaGa_sT_DMSO_EV","WaGa_RNA_147","X2706_WaGa_wt_EV")

    col_order <- c("gene_name",  "MKL.1_RNA","MKL.1_RNA_118","MKL.1_RNA_147","MKL.1_EV.RNA","MKL.1_EV.RNA_2","MKL.1_EV.RNA_118","MKL.1_EV.RNA_87","MKL.1_EV.RNA_27","X042_MKL.1_wt_EV","X042_MKL.1_sT_DMSO","X0505_MKL.1_sT_DMSO_EV","X042_MKL.1_scr_DMSO_EV","X0505_MKL.1_scr_DMSO_EV","X042_MKL.1_sT_Dox","X0505_MKL.1_sT_Dox_EV","X042_MKL.1_scr_Dox_EV","X0505_MKL.1_scr_Dox_EV",     "Geneid.1","gene_name.1",    "WaGa_RNA","WaGa_RNA_118","WaGa_RNA_147",    "WaGa_EV.RNA","WaGa_EV.RNA_2","WaGa_EV.RNA_118","WaGa_EV.RNA_147","WaGa_EV.RNA_226","X1107_WaGa_wt_EV","X1605_WaGa_wt_EV","X2706_WaGa_wt_EV",    "X1107_WaGa_sT_DMSO_EV","X1605_WaGa_sT_DMSO_EV","X2706_WaGa_sT_DMSO_EV",    "X1107_WaGa_scr_DMSO_EV","X1605_WaGa_scr_DMSO_EV","X2706_WaGa_scr_DMSO_EV",     "X1107_WaGa_sT_Dox_EV","X1605_WaGa_sT_Dox_EV","X2706_WaGa_sT_Dox_EV",     "X1107_WaGa_scr_Dox_EV","X1605_WaGa_scr_Dox_EV","X2706_WaGa_scr_Dox_EV")
    reordered.raw_human <- d.raw_human[,col_order]

    d.raw_virus <- read.delim2("../merged_gene_counts_virus_rounded.txt",sep="\t", header=TRUE, row.names=1)
    reordered.raw_virus <- d.raw_virus[,col_order]

    identical(colnames(reordered.raw_human), colnames(reordered.raw_virus))

    reordered.raw <- rbind(reordered.raw_human, reordered.raw_virus)

    #rename
    colnames(reordered.raw) <- c("gene_name", "MKL-1 parental cell RNA","MKL-1 parental cell RNA 118","MKL-1 parental cell RNA 147",    "MKL-1 wt EV RNA","MKL-1 wt EV RNA 2","MKL-1 wt EV RNA 118","MKL-1 wt EV RNA 87","MKL-1 wt EV RNA 27","MKL-1 wt EV RNA 042",    "MKL-1 sT DMSO EV RNA 042","MKL-1 sT DMSO EV RNA 0505",  "MKL-1 scr DMSO EV RNA 042","MKL-1 scr DMSO EV RNA 0505",  "MKL-1 sT Dox EV RNA 042","MKL-1 sT Dox EV RNA 0505",  "MKL-1 scr Dox EV RNA 042","MKL-1 scr Dox EV RNA 0505",   "Geneid.1","gene_name.1",    "WaGa parental cell RNA","WaGa parental cell RNA 118","WaGa parental cell RNA 147",    "WaGa wt EV RNA","WaGa wt EV RNA 2","WaGa wt EV RNA 118","WaGa wt EV RNA 147","WaGa wt EV RNA 226","WaGa wt EV RNA 1107","WaGa wt EV RNA 1605","WaGa wt EV RNA 2706",    "WaGa sT DMSO EV RNA 1107","WaGa sT DMSO EV RNA 1605","WaGa sT DMSO EV RNA 2706",     "WaGa scr DMSO EV RNA 1107","WaGa scr DMSO EV RNA 1605","WaGa scr DMSO EV RNA 2706",     "WaGa sT Dox EV RNA 1107","WaGa sT Dox EV RNA 1605","WaGa sT Dox EV RNA 2706",     "WaGa scr Dox EV RNA 1107","WaGa scr Dox EV RNA 1605","WaGa scr Dox EV RNA 2706")

    reordered.raw$gene_name <- NULL
    reordered.raw$Geneid.1 <- NULL
    reordered.raw$gene_name.1 <- NULL
    write.csv(reordered.raw, file="counts.txt")
    #IMPORTANT that we should filter the data with the counts in the STEP!
    d <- reordered.raw[rowSums(reordered.raw>3)>2,]

    condition_for_pca = as.factor(c("RNA","RNA","RNA","EV","EV","EV","EV","EV","EV","sT.DMSO","sT.DMSO","scr.DMSO","scr.DMSO","sT.Dox","sT.Dox","scr.Dox","scr.Dox",    "RNA","RNA","RNA","EV","EV","EV","EV","EV","EV","EV","EV","sT.DMSO","sT.DMSO","sT.DMSO","scr.DMSO","scr.DMSO","scr.DMSO","sT.Dox","sT.Dox","sT.Dox","scr.Dox","scr.Dox","scr.Dox"))

    condition = as.factor(c("MKL1.RNA","MKL1.RNA","MKL1.RNA","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.sT.DMSO","MKL1.sT.DMSO","MKL1.scr.DMSO","MKL1.scr.DMSO","MKL1.sT.Dox","MKL1.sT.Dox","MKL1.scr.Dox","MKL1.scr.Dox",    "WaGa.RNA","WaGa.RNA","WaGa.RNA","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.sT.DMSO","WaGa.sT.DMSO","WaGa.sT.DMSO","WaGa.scr.DMSO","WaGa.scr.DMSO","WaGa.scr.DMSO","WaGa.sT.Dox","WaGa.sT.Dox","WaGa.sT.Dox","WaGa.scr.Dox","WaGa.scr.Dox","WaGa.scr.Dox"))

    donor = as.factor(c("1","118","147",  "1","2","118","87","27","042",  "042","0505","042","0505","042","0505","042","0505",    "1","118","147",  "1","2","118","147","226","1107","1605","2706",   "1107","1605","2706","1107","1605","2706","1107","1605","2706","1107","1605","2706"))

    batch = as.factor(c("2021.08","2021.09","2021.09","2021.08","2021.08","2021.09","2021.09","2021.09","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08",    "2021.08","2021.09","2021.09","2021.08","2021.08","2021.09","2021.09","2021.09","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11"))

    cell.line = as.factor(c("MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1",    "WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa"))

    ids = as.factor(c("MKL1.RNA","MKL1.RNA.118","MKL1.RNA.147",  "MKL1.EV","MKL1.EV.2","MKL1.EV.118","MKL1.EV.87","MKL1.EV.27","MKL1.EV.042",  "MKL1.EV.sT.DMSO.042","MKL1.EV.sT.DMSO.0505",  "MKL1.EV.scr.DMSO.042","MKL1.EV.scr.DMSO.0505",  "MKL1.EV.sT.Dox.042","MKL1.EV.sT.Dox.0505",  "MKL1.EV.scr.Dox.042","MKL1.EV.scr.Dox.0505",      "WaGa.RNA","WaGa.RNA.118","WaGa.RNA.147",  "WaGa.EV","WaGa.EV.2","WaGa.EV.118","WaGa.EV.147","WaGa.EV.226","WaGa.EV.1107","WaGa.EV.1605","WaGa.EV.2706",  "WaGa.EV.sT.DMSO.1107","WaGa.EV.sT.DMSO.1605","WaGa.EV.sT.DMSO.2706",  "WaGa.EV.scr.DMSO.1107","WaGa.EV.scr.DMSO.1605","WaGa.EV.scr.DMSO.2706",  "WaGa.EV.sT.Dox.1107","WaGa.EV.sT.Dox.1605","WaGa.EV.sT.Dox.2706",  "WaGa.EV.scr.Dox.1107","WaGa.EV.scr.Dox.1605","WaGa.EV.scr.Dox.2706"))

    cData = data.frame(row.names=colnames(d), condition=condition,  donor=donor, batch=batch, cell.line=cell.line, ids=ids)
    dds<-DESeqDataSetFromMatrix(countData=d, colData=cData, design=~batch+condition)

    #rld <- rlogTransformation(dds)
    rld <- vst(dds)

Preparing the data for PCA_MKL1 and PCA_WaGa drawing

    #d_WaGa <- d[, !grepl("parental|MKL-1", names(d))]
    d_WaGa <- d[, !grepl("MKL-1", names(d))]
    condition = as.factor(c("WaGa.RNA","WaGa.RNA","WaGa.RNA","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.EV","WaGa.sT.DMSO","WaGa.sT.DMSO","WaGa.sT.DMSO","WaGa.scr.DMSO","WaGa.scr.DMSO","WaGa.scr.DMSO","WaGa.sT.Dox","WaGa.sT.Dox","WaGa.sT.Dox","WaGa.scr.Dox","WaGa.scr.Dox","WaGa.scr.Dox"))
    donor = as.factor(c("1","118","147","1","2","118","147","226","1107","1605","2706",   "1107","1605","2706","1107","1605","2706","1107","1605","2706","1107","1605","2706"))
    batch = as.factor(c("2021.08","2021.09","2021.09","2021.08","2021.08","2021.09","2021.09","2021.09","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11","2022.11"))
    cell.line = as.factor(c("WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa","WaGa"))
    ids = as.factor(c("WaGa.RNA","WaGa.RNA.118","WaGa.RNA.147","WaGa.EV","WaGa.EV.2","WaGa.EV.118","WaGa.EV.147","WaGa.EV.226","WaGa.EV.1107","WaGa.EV.1605","WaGa.EV.2706",  "WaGa.EV.sT.DMSO.1107","WaGa.EV.sT.DMSO.1605","WaGa.EV.sT.DMSO.2706",  "WaGa.EV.scr.DMSO.1107","WaGa.EV.scr.DMSO.1605","WaGa.EV.scr.DMSO.2706",  "WaGa.EV.sT.Dox.1107","WaGa.EV.sT.Dox.1605","WaGa.EV.sT.Dox.2706",  "WaGa.EV.scr.Dox.1107","WaGa.EV.scr.Dox.1605","WaGa.EV.scr.Dox.2706"))
    cData = data.frame(row.names=colnames(d_WaGa), condition=condition,  donor=donor, batch=batch, cell.line=cell.line, ids=ids)
    dds_WaGa<-DESeqDataSetFromMatrix(countData=d_WaGa, colData=cData, design=~batch+condition)
    rld_WaGa <- vst(dds_WaGa)

    #d_MKL1 <- d[, !grepl("parental|WaGa", names(d))]
    d_MKL1 <- d[, !grepl("WaGa", names(d))]
    condition = as.factor(c("MKL1.RNA","MKL1.RNA","MKL1.RNA","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.EV","MKL1.sT.DMSO","MKL1.sT.DMSO","MKL1.scr.DMSO","MKL1.scr.DMSO","MKL1.sT.Dox","MKL1.sT.Dox","MKL1.scr.Dox","MKL1.scr.Dox"))
    donor = as.factor(c("1","118","147","1","2","118","87","27","042",  "042","0505","042","0505","042","0505","042","0505"))
    batch = as.factor(c("2021.08","2021.09","2021.09","2021.08","2021.08","2021.09","2021.09","2021.09","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08","2022.08"))
    cell.line = as.factor(c("MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1","MKL-1"))
    ids = as.factor(c("MKL1.RNA","MKL1.RNA.118","MKL1.RNA.147", "MKL1.EV","MKL1.EV.2","MKL1.EV.118","MKL1.EV.87","MKL1.EV.27","MKL1.EV.042",  "MKL1.EV.sT.DMSO.042","MKL1.EV.sT.DMSO.0505",  "MKL1.EV.scr.DMSO.042","MKL1.EV.scr.DMSO.0505",  "MKL1.EV.sT.Dox.042","MKL1.EV.sT.Dox.0505",  "MKL1.EV.scr.Dox.042","MKL1.EV.scr.Dox.0505"))
    cData = data.frame(row.names=colnames(d_MKL1), condition=condition,  donor=donor, batch=batch, cell.line=cell.line, ids=ids)
    dds_MKL1<-DESeqDataSetFromMatrix(countData=d_MKL1, colData=cData, design=~batch+condition)
    rld_MKL1 <- vst(dds_MKL1)

    # -- before pca --
    png("pca.png", 1200, 800)
    plotPCA(rld, intgroup=c("condition"))
    #plotPCA(rld, intgroup = c("condition", "batch"))
    #plotPCA(rld, intgroup = c("condition", "ids"))
    #plotPCA(rld, "batch")
    dev.off()
    #72% (PC1), 11% (PC2)   8% (PC3)

    png("pca_WaGa.png", 1200, 800)
    plotPCA(rld_WaGa, intgroup = c("condition"))
    dev.off()
    #80%, 9%

    png("pca_MKL1.png", 1200, 800)
    plotPCA(rld_MKL1, intgroup = c("condition"))
    dev.off()
    #77%, 15%

PCA plot untreated/wt vs parental cells; 1x für WaGa cell line und 1x für MKL-1 cells

    #install.packages("BiocManager")
    #BiocManager::install("genefilter")
    #install.packages("writexl")
    #install.packages("plotly")
    #install.packages("dplyr")
    #install.packages("ggrepel")
    # Load required libraries
    library(ggrepel)
    library(dplyr)
    library(ggplot2)
    library(genefilter)
    library(writexl)

    # ---- drawing PCA WaGa ----

    # Perform PCA on WaGa dataset
    pca_data <- plotPCA(rld_WaGa, intgroup = c("condition", "batch", "cell.line"), returnData = TRUE)
    write.csv(pca_data, file = "PCA_data_WaGa.csv")

    # Compute top variable genes and perform PCA
    top_variable_genes <- 500
    gene_variance <- rowVars(assay(rld_WaGa))
    selected_genes <- order(gene_variance, decreasing = TRUE)[seq_len(min(top_variable_genes, length(gene_variance)))]

    # Extract data and run PCA
    expression_matrix <- t(assay(rld_WaGa)[selected_genes, ])
    pca_results <- prcomp(expression_matrix, center = TRUE, scale. = TRUE)
    summary(pca_results)  # Display variance explained by PCs

    # Extract variance percentages
    variance_explained <- summary(pca_results)$importance[2, ]
    pc1_variance <- sprintf("%.1f", variance_explained[1] * 100)
    pc2_variance <- sprintf("%.1f", variance_explained[2] * 100)

    # Create a PCA dataframe
    pca_df <- as.data.frame(pca_results$x)
    rownames(pca_df) <- rownames(pca_data)

    # Merge PCA results with original metadata
    pca_data <- pca_data %>%
    select(-PC1, -PC2)  # Remove old PCA values

    merged_pca_data <- merge(pca_data, pca_df, by = "row.names")
    rownames(merged_pca_data) <- merged_pca_data$Row.names
    merged_pca_data$Row.names <- NULL

    # Save merged PCA data
    write.csv(merged_pca_data, file = "Merged_PCA_WaGa.csv")
    write_xlsx(merged_pca_data, path = "Merged_PCA_WaGa.xlsx")

    # Define condition labels
    condition_labels <- c(
          "WaGa.EV" = "wt EV RNA",
          "WaGa.scr.DMSO" = "scr DMSO EV RNA",
          "WaGa.scr.Dox" = "scr Dox EV RNA",
          "WaGa.sT.DMSO" = "sT DMSO EV RNA",
          "WaGa.sT.Dox" = "sT Dox EV RNA",
          "WaGa.RNA" = "parental cell RNA"
    )

    # Apply condition labels and filter relevant conditions
    filtered_pca_data <- merged_pca_data %>%
    mutate(condition = recode(condition, !!!condition_labels)) %>%
    filter(condition %in% c("wt EV RNA", "parental cell RNA"))

    # Prepare PCA data for plotting
    plot_pca_df <- as.data.frame(filtered_pca_data[, c("PC1", "PC2")])
    plot_pca_df$condition <- filtered_pca_data$condition
    plot_pca_df$sample_name <- rownames(filtered_pca_data)

    # Save PCA plot as PNG
    png("PCA_WaGa.png", width = 1000, height = 600, res = 150)

    ggplot(plot_pca_df, aes(x = PC1, y = PC2, color = condition)) +
    geom_point(size = 4, alpha = 0.8) +  # Plot points
    labs(
    title = "",
    x = paste0("Principal Component 1 (", pc1_variance, "%)"),
    y = paste0("Principal Component 2 (", pc2_variance, "%)")
    ) +
    theme_minimal() +
    scale_color_manual(values = c("blue", "red"))  # Customize colors

    dev.off()  # Close the PNG device and save the file

    # ---- drawing PCA MKL-1 ----

    # Compute PCA on MKL1 dataset
    pca_data <- plotPCA(rld_MKL1, intgroup = c("condition", "batch", "cell.line"), returnData = TRUE)
    write.csv(pca_data, file = "PCA_data_MKL1.csv")

    # Compute top variable genes and perform PCA
    top_variable_genes <- 500
    gene_variance <- rowVars(assay(rld_MKL1))
    selected_genes <- order(gene_variance, decreasing = TRUE)[seq_len(min(top_variable_genes, length(gene_variance)))]

    # Extract and transpose data for PCA
    expression_matrix <- t(assay(rld_MKL1)[selected_genes, ])
    pca_results <- prcomp(expression_matrix, center = TRUE, scale. = TRUE)

    # Extract variance explained percentages
    variance_explained <- summary(pca_results)$importance[2, ]
    pc1_variance <- sprintf("%.1f", variance_explained[1] * 100)
    pc2_variance <- sprintf("%.1f", variance_explained[2] * 100)

    # Create a PCA dataframe
    pca_df <- as.data.frame(pca_results$x)
    rownames(pca_df) <- rownames(pca_data)

    # Remove previous PC1 & PC2 columns and merge new PCA results
    pca_data <- pca_data %>%
    select(-PC1, -PC2)

    merged_pca_data <- merge(pca_data, pca_df, by = "row.names")
    rownames(merged_pca_data) <- merged_pca_data$Row.names
    merged_pca_data$Row.names <- NULL

    # Save merged PCA data
    write.csv(merged_pca_data, file = "Merged_PCA_MKL1.csv")
    write_xlsx(merged_pca_data, path = "Merged_PCA_MKL1.xlsx")

    # Define condition labels for clarity
    condition_labels <- c(
    "MKL1.EV" = "wt EV RNA",
    "MKL1.scr.DMSO" = "scr DMSO EV RNA",
    "MKL1.scr.Dox" = "scr Dox EV RNA",
    "MKL1.sT.DMSO" = "sT DMSO EV RNA",
    "MKL1.sT.Dox" = "sT Dox EV RNA",
    "MKL1.RNA" = "parental cell RNA"
    )

    # Apply condition labels and filter for relevant conditions
    filtered_pca_data <- merged_pca_data %>%
    mutate(condition = recode(condition, !!!condition_labels)) %>%
    filter(condition %in% c("wt EV RNA", "parental cell RNA"))

    # Prepare PCA data for plotting
    plot_pca_df <- as.data.frame(filtered_pca_data[, c("PC1", "PC2")])
    plot_pca_df$condition <- filtered_pca_data$condition
    plot_pca_df$sample_name <- rownames(filtered_pca_data)

    # Save PCA plot as PNG
    png("PCA_MKL1.png", width = 1000, height = 600, res = 150)

    ggplot(plot_pca_df, aes(x = PC1, y = PC2, color = condition)) +
    geom_point(size = 4, alpha = 0.8) +  # Plot points
    labs(
    title = "",
    x = paste0("Principal Component 1 (", pc1_variance, "%)"),
    y = paste0("Principal Component 2 (", pc2_variance, "%)")
    ) +
    theme_minimal() +
    scale_color_manual(values = c("blue", "red"))  # Customize colors

    dev.off()  # Close the PNG device and save the file

Convert bam to bigwig using deepTools by feeding inverse of DESeq’s size Factor

    sizeFactors(dds)
    #NULL
    dds <- estimateSizeFactors(dds)
    > sizeFactors(dds)
    #WaGa parental cell RNA*
    #                  1.1242096                   2.1097851
    #WaGa parental cell RNA 118*  WaGa parental cell RNA 147*
    #                  1.6925780                   1.6712182

    raw_counts <- counts(dds)
    normalized_counts <- counts(dds, normalized=TRUE)
    write.table(raw_counts, file="raw_counts.txt", sep="\t", quote=F, col.names=NA)
    write.table(normalized_counts, file="normalized_counts.txt", sep="\t", quote=F, col.names=NA)

    #1/2,1097851=0,473981924 (!=0.4958732)
    bamCoverage --bam ../markDuplicates/WaGa_RNAAligned.sortedByCoord.out.markDups.bam -o WaGa_RNA.bw --binSize 10 --scaleFactor 0.473981924 --effectiveGenomeSize 2864785220
    #1/1,6925780=0,590814722 (!=0.6013898)
    bamCoverage --bam ../markDuplicates/WaGa_RNA_118Aligned.sortedByCoord.out.markDups.bam -o WaGa_RNA_118.bw --binSize 10 --scaleFactor 0.590814722 --effectiveGenomeSize 2864785220
    #1/1,6712182=0,598365911 (!=0.6154516)
    bamCoverage --bam ../markDuplicates/WaGa_RNA_147Aligned.sortedByCoord.out.markDups.bam -o WaGa_RNA_147.bw --binSize 10 --scaleFactor 0.598365911 --effectiveGenomeSize 2864785220

DEG analysis

    ##https://bioconductor.statistik.tu-dortmund.de/packages/3.7/data/annotation/
    #BiocManager::install("EnsDb.Mmusculus.v79")
    #library(EnsDb.Mmusculus.v79)
    #edb <- EnsDb.Mmusculus.v79
    #https://bioconductor.org/packages/release/bioc/vignettes/biomaRt/inst/doc/accessing_ensembl.html#selecting-an-ensembl-biomart-database-and-dataset
    #https://bioconductor.org/packages/release/bioc/vignettes/biomaRt/inst/doc/accessing_ensembl.html#selecting-an-ensembl-biomart-database-and-dataset
    library(biomaRt)
    #listMarts() == listEnsembl() == Ensembl Genes 113
    listEnsembl()
    listMarts()
    listEnsemblArchives()

    #mkdir degenes
    setwd("degenes")

    dds$condition <- relevel(dds$condition, "WaGa.RNA")
    dds = DESeq(dds, betaPrior=FALSE)
    resultsNames(dds)
    clist <- c("WaGa.EV_vs_WaGa.RNA")

    dds$condition <- relevel(dds$condition, "MKL1.RNA")
    dds = DESeq(dds, betaPrior=FALSE)
    resultsNames(dds)
    clist <- c("MKL1.EV_vs_MKL1.RNA")

    #ensembl <- useEnsembl(biomart = "ensembl", dataset = "hsapiens_gene_ensembl", version="112")
    #ensembl <- useMart("ensembl", dataset="hsapiens_gene_ensembl", host="https://may2023.archive.ensembl.org")
    ensembl <- useMart("ensembl", dataset="hsapiens_gene_ensembl", host="https://www.ensembl.org")

    listDatasets(ensembl)
    attributes = listAttributes(ensembl)
    attributes[1:25,]

    for (i in clist) {
          #i<-clist[1]
          contrast = paste("condition", i, sep="_")
          res = results(dds, name=contrast)
          res <- res[!is.na(res$log2FoldChange),]
          geness <- getBM(attributes = c('ensembl_gene_id', 'external_gene_name', 'gene_biotype', 'entrezgene_id', 'chromosome_name', 'start_position', 'end_position', 'strand', 'description'), filters = 'ensembl_gene_id', values = rownames(res), mart = ensembl)
          geness_uniq <- distinct(geness, ensembl_gene_id, .keep_all= TRUE)

          #merge by column by common colunmn name, in the case "GENEID"
          res$ENSEMBL = rownames(res)
          identical(rownames(res), rownames(geness_uniq))
          res_df <- as.data.frame(res)
          geness_res <- merge(geness_uniq, res_df, by.x="ensembl_gene_id", by.y="ENSEMBL")
          dim(geness_res)
          rownames(geness_res) <- geness_res$ensembl_gene_id
          geness_res$ensembl_gene_id <- NULL
          write.csv(as.data.frame(geness_res[order(geness_res$pvalue),]), file = paste(i, "all.txt", sep="-"))
          up <- subset(geness_res, padj<=0.05 & log2FoldChange>=2)
          down <- subset(geness_res, padj<=0.05 & log2FoldChange<=-2)
          write.csv(as.data.frame(up[order(up$log2FoldChange,decreasing=TRUE),]), file = paste(i, "up.txt", sep="-"))
          write.csv(as.data.frame(down[order(abs(down$log2FoldChange),decreasing=TRUE),]), file = paste(i, "down.txt", sep="-"))
    }

Volcano plot untreated/wt vs parental; 1x für WaGa cell line und 1x für MKL-1 cells

    #A canonical visualization for interpreting differential gene expression results is the volcano plot.
    library(ggrepel)

    # -- WaGa.EV_vs_WaGa.RNA --
    geness_res <- read.csv(file = paste("WaGa.EV_vs_WaGa.RNA", "all.txt", sep="-"), row.names=1)
    geness_res$Color <- "NS or log2FC < 2.0"
    geness_res$Color[geness_res$pvalue < 0.05] <- "P < 0.05"
    geness_res$Color[geness_res$padj < 0.05] <- "P-adj < 0.05"
    geness_res$Color[geness_res$padj < 0.001] <- "P-adj < 0.001"
    geness_res$Color[abs(geness_res$log2FoldChange) < 2.0] <- "NS or log2FC < 2.0"
    geness_res$Color <- factor(geness_res$Color,
                            levels = c("NS or log2FC < 2.0", "P < 0.05",
                                  "P-adj < 0.05", "P-adj < 0.001"))

    geness_res$invert_P <- (-log10(geness_res$pvalue)) * sign(geness_res$log2FoldChange)
    top_g <- c()
    top_g <- c(top_g, geness_res[, 'external_gene_name'][order(geness_res[, 'invert_P'], decreasing = TRUE)[1:200]], geness_res[, 'external_gene_name'][order(geness_res[, 'invert_P'], decreasing = FALSE)[1:200]])
    top_g <- unique(top_g)
    geness_res <- geness_res[, -1*ncol(geness_res)]  #remove invert_P from matrix

    png("WaGa_wt.EV_vs_parental.png",width=1400, height=1000)
    ggplot(geness_res,
          aes(x = log2FoldChange, y = -log10(pvalue),
          color = Color, label = external_gene_name)) +
    geom_vline(xintercept = c(2.0, -2.0), lty = "dashed") +
    geom_hline(yintercept = -log10(0.05), lty = "dashed") +
    geom_point() +
    labs(x = "log2(FC)",
          y = "Significance, -log10(P)",
          color = "Significance") +
    scale_color_manual(values = c(`P-adj < 0.001` = "dodgerblue",
                                  `P-adj < 0.05` = "lightblue",
                                  `P < 0.05` = "orange2",
                                  `NS or log2FC < 2.0` = "gray"),
                      guide = guide_legend(override.aes = list(size = 4))) +
    scale_y_continuous(expand = expansion(mult = c(0,0.05))) +
    geom_text_repel(data = subset(geness_res, external_gene_name %in% top_g & pvalue < 0.05 & (abs(geness_res$log2FoldChange) >= 2.0)),
                      size = 4, point.padding = 0.15, color = "black",
                      min.segment.length = .1, box.padding = .2, lwd = 2) +
    theme_bw(base_size = 16) +
    theme(legend.position = "bottom")
    dev.off()

    # -- MKL1.EV_vs_MKL1.RNA --
    geness_res <- read.csv(file = paste("MKL1.EV_vs_MKL1.RNA", "all.txt", sep="-"), row.names=1)
    geness_res$Color <- "NS or log2FC < 2.0"
    geness_res$Color[geness_res$pvalue < 0.05] <- "P < 0.05"
    geness_res$Color[geness_res$padj < 0.05] <- "P-adj < 0.05"
    geness_res$Color[geness_res$padj < 0.001] <- "P-adj < 0.001"
    geness_res$Color[abs(geness_res$log2FoldChange) < 2.0] <- "NS or log2FC < 2.0"
    geness_res$Color <- factor(geness_res$Color,
                            levels = c("NS or log2FC < 2.0", "P < 0.05",
                                  "P-adj < 0.05", "P-adj < 0.001"))

    geness_res$invert_P <- (-log10(geness_res$pvalue)) * sign(geness_res$log2FoldChange)
    top_g <- c()
    top_g <- c(top_g, geness_res[, 'external_gene_name'][order(geness_res[, 'invert_P'], decreasing = TRUE)[1:200]], geness_res[, 'external_gene_name'][order(geness_res[, 'invert_P'], decreasing = FALSE)[1:200]])
    top_g <- unique(top_g)
    geness_res <- geness_res[, -1*ncol(geness_res)]  #remove invert_P from matrix

    png("MKL-1_wt.EV_vs_parental.png",width=1400, height=1000)
    ggplot(geness_res,
          aes(x = log2FoldChange, y = -log10(pvalue),
          color = Color, label = external_gene_name)) +
    geom_vline(xintercept = c(2.0, -2.0), lty = "dashed") +
    geom_hline(yintercept = -log10(0.05), lty = "dashed") +
    geom_point() +
    labs(x = "log2(FC)",
          y = "Significance, -log10(P)",
          color = "Significance") +
    scale_color_manual(values = c(`P-adj < 0.001` = "dodgerblue",
                                  `P-adj < 0.05` = "lightblue",
                                  `P < 0.05` = "orange2",
                                  `NS or log2FC < 2.0` = "gray"),
                      guide = guide_legend(override.aes = list(size = 4))) +
    scale_y_continuous(expand = expansion(mult = c(0,0.05))) +
    geom_text_repel(data = subset(geness_res, external_gene_name %in% top_g & pvalue < 0.05 & (abs(geness_res$log2FoldChange) >= 2.0)),
                      size = 4, point.padding = 0.15, color = "black",
                      min.segment.length = .1, box.padding = .2, lwd = 2) +
    theme_bw(base_size = 16) +
    theme(legend.position = "bottom")
    dev.off()

    # -- under console --
    mv WaGa.EV_vs_WaGa.RNA-all.txt WaGa_wt.EV_vs_parental-all.txt
    mv WaGa.EV_vs_WaGa.RNA-up.txt WaGa_wt.EV_vs_parental-up.txt
    mv WaGa.EV_vs_WaGa.RNA-down.txt WaGa_wt.EV_vs_parental-down.txt
    mv MKL1.EV_vs_MKL1.RNA-all.txt MKL-1_wt.EV_vs_parental-all.txt
    mv MKL1.EV_vs_MKL1.RNA-up.txt MKL-1_wt.EV_vs_parental-up.txt
    mv MKL1.EV_vs_MKL1.RNA-down.txt MKL-1_wt.EV_vs_parental-down.txt
    for cmp in "WaGa_wt.EV_vs_parental" "MKL-1_wt.EV_vs_parental"; do
          echo "~/Tools/csv2xls-0.4/csv_to_xls.py ${cmp}-all.txt ${cmp}-up.txt ${cmp}-down.txt -d$',' -o ${cmp}.xls"
    done
    ~/Tools/csv2xls-0.4/csv_to_xls.py WaGa_wt.EV_vs_parental-all.txt WaGa_wt.EV_vs_parental-up.txt WaGa_wt.EV_vs_parental-down.txt -d$',' -o WaGa_wt.EV_vs_parental.xls
    ~/Tools/csv2xls-0.4/csv_to_xls.py MKL-1_wt.EV_vs_parental-all.txt MKL-1_wt.EV_vs_parental-up.txt MKL-1_wt.EV_vs_parental-down.txt -d$',' -o MKL-1_wt.EV_vs_parental.xls

Heatmap untreated/wt vs parental; 1x für WaGa cell line und 1x für MKL-1 cells

    install.packages("gplots")
    library("gplots")

    # -- WaGa cell line --
    cut -d',' -f1-1 ./WaGa_wt.EV_vs_parental-up.txt > WaGa_wt.EV_vs_parental-up.id
    cut -d',' -f1-1 ./WaGa_wt.EV_vs_parental-down.txt > WaGa_wt.EV_vs_parental-down.id
    cat WaGa_wt.EV_vs_parental-up.id WaGa_wt.EV_vs_parental-down.id | sort -u > ids  #24207

    #add Gene_Id in the first line.
    GOI <- read.csv("ids")$Gene_Id
    RNASeq.NoCellLine <- assay(rld)

    #clustering methods: "ward.D", "ward.D2", "single", "complete", "average" (= UPGMA), "mcquitty" (= WPGMA), "median" (= WPGMC) or "centroid" (= UPGMC).  pearson or spearman
    datamat = RNASeq.NoCellLine[GOI, c("WaGa parental cell RNA","WaGa parental cell RNA 118","WaGa parental cell RNA 147", "WaGa wt EV RNA","WaGa wt EV RNA 2","WaGa wt EV RNA 118","WaGa wt EV RNA 147","WaGa wt EV RNA 226",  "WaGa wt EV RNA 1107","WaGa wt EV RNA 1605","WaGa wt EV RNA 2706")]
    hr <- hclust(as.dist(1-cor(t(datamat), method="pearson")), method="average")
    hc <- hclust(as.dist(1-cor(datamat, method="spearman")), method="average")
    mycl = cutree(hr, h=max(hr$height)/1.05)
    mycol = c("YELLOW", "DARKBLUE", "DARKORANGE", "DARKMAGENTA", "DARKCYAN", "DARKRED",  "MAROON", "DARKGREEN", "LIGHTBLUE", "PINK", "MAGENTA", "LIGHTCYAN","LIGHTGREEN", "BLUE", "ORANGE", "CYAN", "RED", "GREEN");

    mycol = mycol[as.vector(mycl)]
    sampleCols <- rep('GREY',ncol(datamat))
    names(sampleCols) <- c("WaGa parental cell RNA","WaGa parental cell RNA 118","WaGa parental cell RNA 147", "WaGa wt EV RNA","WaGa wt EV RNA 2","WaGa wt EV RNA 118","WaGa wt EV RNA 147","WaGa wt EV RNA 226",  "WaGa wt EV RNA 1107","WaGa wt EV RNA 1605","WaGa wt EV RNA 2706")

    #DARKBLUE DARKGREEN DARKRED DARKORANGE
    sampleCols["WaGa parental cell RNA"] <- 'GREY'
    sampleCols["WaGa parental cell RNA 118"] <- 'GREY'
    sampleCols["WaGa parental cell RNA 147"] <- 'GREY'
    sampleCols["WaGa wt EV RNA"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 2"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 118"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 147"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 226"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 1107"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 1605"] <- 'GREEN'
    sampleCols["WaGa wt EV RNA 2706"] <- 'GREEN'

    png("DEGs_Heatmap_WaGa.png", width=1000, height=1200)
    heatmap.2(as.matrix(datamat),Rowv=as.dendrogram(hr),Colv = NA, dendrogram = 'row',
                scale='row',trace='none',col=bluered(75),
                RowSideColors = mycol, ColSideColors = sampleCols, labRow="", margins=c(22,10), cexRow=8, cexCol=2, srtCol=30, lwid=c(1,7), lhei = c(1, 8), legend("top", title = "",legend=c("WaGa parental cell RNA","WaGa wt EV RNA"), fill=c("GREY","GREEN"), cex=1.0, box.lty=0))
    dev.off()

    # save cluster members
    subset_1<-names(subset(mycl, mycl == '1'))
    subset_1_ <- getBM(attributes = c('ensembl_gene_id', 'external_gene_name', 'gene_biotype', 'entrezgene_id', 'chromosome_name', 'start_position', 'end_position', 'strand', 'description'),
          filters = 'ensembl_gene_id',
          values = subset_1,
          mart = ensembl)
    subset_1_uniq <- distinct(subset_1_, ensembl_gene_id, .keep_all= TRUE)
    subset_1_expr  <- datamat[subset_1,]
    subset_1_expr <- as.data.frame(subset_1_expr)
    subset_1_expr$ENSEMBL = rownames(subset_1_expr)
    cluster1_YELLOW <- merge(subset_1_uniq, subset_1_expr, by.x="ensembl_gene_id", by.y="ENSEMBL")
    #write.csv(cluster1_YELLOW,file='cluster1_YELLOW.txt')

    subset_2<-names(subset(mycl, mycl == '2'))
    subset_2_ <- getBM(attributes = c('ensembl_gene_id', 'external_gene_name', 'gene_biotype', 'entrezgene_id', 'chromosome_name', 'start_position', 'end_position', 'strand', 'description'),
          filters = 'ensembl_gene_id',
          values = subset_2,
          mart = ensembl)
    subset_2_uniq <- distinct(subset_2_, ensembl_gene_id, .keep_all= TRUE)
    subset_2_expr  <- datamat[subset_2,]
    subset_2_expr <- as.data.frame(subset_2_expr)
    subset_2_expr$ENSEMBL = rownames(subset_2_expr)
    cluster2_DARKBLUE <- merge(subset_2_uniq, subset_2_expr, by.x="ensembl_gene_id", by.y="ENSEMBL")
    #write.csv(cluster2_DARKBLUE,file='cluster2_DARKBLUE.txt')

    write_xlsx(list(
    "Cluster 1 YELLOW" = cluster1_YELLOW,
    "Cluster 2 DARKBLUE" = cluster2_DARKBLUE
    ), "DEGs_heatmap_data_WaGa.xlsx")

    # -- MKL-1 cell line --
    cut -d',' -f1-1 ./MKL-1_wt.EV_vs_parental-up.txt > MKL-1_wt.EV_vs_parental-up.id
    cut -d',' -f1-1 ./MKL-1_wt.EV_vs_parental-down.txt > MKL-1_wt.EV_vs_parental-down.id
    cat MKL-1_wt.EV_vs_parental-up.id MKL-1_wt.EV_vs_parental-down.id | sort -u > ids  #20720

    #add Gene_Id in the first line.
    GOI <- read.csv("ids")$Gene_Id
    RNASeq.NoCellLine <- assay(rld)

    #clustering methods: "ward.D", "ward.D2", "single", "complete", "average" (= UPGMA), "mcquitty" (= WPGMA), "median" (= WPGMC) or "centroid" (= UPGMC).  pearson or spearman
    datamat = RNASeq.NoCellLine[GOI, c("MKL-1 parental cell RNA","MKL-1 parental cell RNA 118","MKL-1 parental cell RNA 147", "MKL-1 wt EV RNA","MKL-1 wt EV RNA 2","MKL-1 wt EV RNA 27","MKL-1 wt EV RNA 87","MKL-1 wt EV RNA 118",  "MKL-1 wt EV RNA 042")]
    #Check for Zero Variance Rows; If any row has variance 0, remove it.
    #datamat <- datamat[apply(datamat, 1, var) != 0, ]
    hr <- hclust(as.dist(1-cor(t(datamat), method="pearson")), method="average")
    hc <- hclust(as.dist(1-cor(datamat, method="spearman")), method="average")
    mycl = cutree(hr, h=max(hr$height)/1.05)
    mycol = c("YELLOW", "DARKBLUE", "DARKORANGE", "DARKMAGENTA", "DARKCYAN", "DARKRED",  "MAROON", "DARKGREEN", "LIGHTBLUE", "PINK", "MAGENTA", "LIGHTCYAN","LIGHTGREEN", "BLUE", "ORANGE", "CYAN", "RED", "GREEN");

    mycol = mycol[as.vector(mycl)]
    sampleCols <- rep('GREY',ncol(datamat))
    names(sampleCols) <- c("MKL-1 parental cell RNA","MKL-1 parental cell RNA 118","MKL-1 parental cell RNA 147", "MKL-1 wt EV RNA","MKL-1 wt EV RNA 2","MKL-1 wt EV RNA 27","MKL-1 wt EV RNA 87","MKL-1 wt EV RNA 118",  "MKL-1 wt EV RNA 042")

    sampleCols["MKL-1 parental cell RNA"] <- 'GREY'
    sampleCols["MKL-1 parental cell RNA 118"] <- 'GREY'
    sampleCols["MKL-1 parental cell RNA 147"] <- 'GREY'
    sampleCols["MKL-1 wt EV RNA"] <- 'GREEN'
    sampleCols["MKL-1 wt EV RNA 2"] <- 'GREEN'
    sampleCols["MKL-1 wt EV RNA 27"] <- 'GREEN'
    sampleCols["MKL-1 wt EV RNA 87"] <- 'GREEN'
    sampleCols["MKL-1 wt EV RNA 118"] <- 'GREEN'
    sampleCols["MKL-1 wt EV RNA 042"] <- 'GREEN'

    png("DEGs_Heatmap_MKL-1.png", width=1000, height=1200)
    heatmap.2(as.matrix(datamat),Rowv=as.dendrogram(hr),Colv = NA, dendrogram = 'row',
                scale='row',trace='none',col=bluered(75),
                RowSideColors = mycol, ColSideColors = sampleCols, labRow="", margins=c(22,10), cexRow=8, cexCol=2, srtCol=30, lwid=c(1,7), lhei = c(1, 8), legend("top", title = "",legend=c("MKL-1 parental cell RNA","MKL-1 wt EV RNA"), fill=c("GREY","GREEN"), cex=1.0, box.lty=0))
    dev.off()

    # save cluster members
    subset_1<-names(subset(mycl, mycl == '1'))
    subset_1_ <- getBM(attributes = c('ensembl_gene_id', 'external_gene_name', 'gene_biotype', 'entrezgene_id', 'chromosome_name', 'start_position', 'end_position', 'strand', 'description'),
          filters = 'ensembl_gene_id',
          values = subset_1,
          mart = ensembl)
    subset_1_uniq <- distinct(subset_1_, ensembl_gene_id, .keep_all= TRUE)
    subset_1_expr  <- datamat[subset_1,]
    subset_1_expr <- as.data.frame(subset_1_expr)
    subset_1_expr$ENSEMBL = rownames(subset_1_expr)
    cluster1_YELLOW <- merge(subset_1_uniq, subset_1_expr, by.x="ensembl_gene_id", by.y="ENSEMBL")
    #write.csv(cluster1_YELLOW,file='cluster1_YELLOW.txt')

    subset_2<-names(subset(mycl, mycl == '2'))
    subset_2_ <- getBM(attributes = c('ensembl_gene_id', 'external_gene_name', 'gene_biotype', 'entrezgene_id', 'chromosome_name', 'start_position', 'end_position', 'strand', 'description'),
          filters = 'ensembl_gene_id',
          values = subset_2,
          mart = ensembl)
    subset_2_uniq <- distinct(subset_2_, ensembl_gene_id, .keep_all= TRUE)
    subset_2_expr  <- datamat[subset_2,]
    subset_2_expr <- as.data.frame(subset_2_expr)
    subset_2_expr$ENSEMBL = rownames(subset_2_expr)
    cluster2_DARKBLUE <- merge(subset_2_uniq, subset_2_expr, by.x="ensembl_gene_id", by.y="ENSEMBL")
    #write.csv(cluster2_DARKBLUE,file='cluster2_DARKBLUE.txt')

    write_xlsx(list(
    "Cluster 1 YELLOW" = cluster1_YELLOW,
    "Cluster 2 DARKBLUE" = cluster2_DARKBLUE
    ), "DEGs_heatmap_data_MKL-1.xlsx")

Distribution of different RNA Species untreated/wt and parental; 1x für WaGa cell line und 1x für MKL-1 cells

    cd ~/DATA/Data_Ute/Data_RNA-Seq_MKL-1+WaGa/results_2025_1/read_distributions
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_RNAAligned.sortedByCoord.out.read_distribution.txt WaGa_parental_cell_RNA_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_RNA_118Aligned.sortedByCoord.out.read_distribution.txt WaGa_parental_cell_RNA_118_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_RNA_147Aligned.sortedByCoord.out.read_distribution.txt WaGa_parental_cell_RNA_147_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_EV-RNAAligned.sortedByCoord.out.read_distribution.txt  WaGa_wt_EV_RNA_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_EV-RNA_2Aligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_2_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_EV-RNA_118Aligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_118_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_EV-RNA_147Aligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_147_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/WaGa_EV-RNA_226Aligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_226_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/1107_WaGa_wt_EVAligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_1107_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/1605_WaGa_wt_EVAligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_1605_read_distribution.txt
    cp ../../../Data_RNA-Seq_WaGa/results/rseqc/read_distribution/2706_WaGa_wt_EVAligned.sortedByCoord.out.read_distribution.txt WaGa_wt_EV_RNA_2706_read_distribution.txt

    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_RNAAligned.sortedByCoord.out.read_distribution.txt MKL-1_parental_cell_RNA_read_distribution.txt
    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_RNA_118Aligned.sortedByCoord.out.read_distribution.txt MKL-1_parental_cell_RNA_118_read_distribution.txt
    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_RNA_147Aligned.sortedByCoord.out.read_distribution.txt MKL-1_parental_cell_RNA_147_read_distribution.txt
    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_EV-RNAAligned.sortedByCoord.out.read_distribution.txt MKL-1_wt_EV_RNA_read_distribution.txt
    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_EV-RNA_2Aligned.sortedByCoord.out.read_distribution.txt MKL-1_wt_EV_RNA_2_read_distribution.txt
    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_EV-RNA_27Aligned.sortedByCoord.out.read_distribution.txt MKL-1_wt_EV_RNA_27_read_distribution.txt
    cp ../../../Data_RNA-Seq_MKL-1/results/rseqc/read_distribution/MKL-1_EV-RNA_87Aligned.sortedByCoord.out.read_distribution.txt MKL-1_wt_EV_RNA_87_read_distribution.txt

RNA fragmentation patterns: is EV-RNA full length or fragmented?

    Since our data consists of single-end RNA-Seq reads, determining the exact RNA fragment length is challenging, as only the read length is available. Unlike paired-end sequencing, where insert size can be calculated, single-end sequencing does not provide direct information about the full fragment size.

Microbial bioinformatics

Microbial bioinformatics uses computational tools to analyze genomes, track evolution, and study functions in microorganisms, including bacteria and viruses.

Leave a Reply Cancel reply