Pipeline for Comparative transcriptomics analysis

Libraries re-mapped

Virus	Project	Libraries	SeqType	Tissue	Ecotype	Genotypes	dpi	Paper	Reference
CaMV	PRJEB49403	6	Paired-End	Rosette	Col-0	WT aphid Mp	14dold, 21dpi	10.1128/spectrum.00136-22	[1]
TuYV	PRJEB49403	3	Paired-End	Rosette	Col-0	WT aphid Mp	14dold, 21dpi	10.1128/spectrum.00136-22	[1]
TuMV	PRJNA788379	8	Paired-End	Rosette	Col-0 and Bar-0	Wt	28dold, 14dpi	10.1371/journal.pone.0275588	[2]
TVCV	PRJNA788379	8	Paired-End	Rosette	Col-0 and Bar-0	Wt	28dold, 14dpi	10.1371/journal.pone.0275588	[2]
TCV	PRJNA336058	6 (12)	Paired-End	Leaves	Col-	WT and dcl1-/-	28dold, 8dpi	10.1038/srep36007.	[3]
CMV	PRJNA1124548	6	Paired-End	Rosette	Col-0	WT, CMV and CMV-Δ2b	14dpi	10.1038/s41467-025-65355-1	[4]
CalCuV	PRJNA244971	8 (20)	Paired-End	Leaves	Col-0	WT	7 leaf stage, 10 and 20dpi	10.1038/nature14171	[5]

Code Running

From FeatureCouns

 Rscript Codes/deseqDEG.R FeatureCounts/Liu_CMV_fcounts.txt "WTvsWT" "CMVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Gyula_TVCV_TuMV_fcounts.txt "WTvsWT" "TuMVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Gyula_TVCV_TuMV_fcounts_mTVCVMockR1.txt "WTvsWT" "TVCVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Gao_TCV_fcounts.txt  "WTvsWT" "TCVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Chesnais_CaMV_TuYV_fcounts.txt  "WTvsWT" "CaMVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Chesnais_CaMV_TuYV_fcounts.txt  "WTvsWT" "TuYVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Zorzatto_CaLCuV_fcounts2.txt "20dpivs20dpi" "CaLCuVvsMock"
 Rscript Codes/deseqDEG.R FeatureCounts/Zorzatto_CaLCuV_fcounts2.txt "10dpivs10dpi" "CaLCuVvsMock"

Gene selection

(cat Data/deseq_head.txt && for gene in `cat Data/atg_genes_DH.txt | cut -f 1` ; do grep $gene FeatureCounts/*DESeq2*; done  )  > Data/ATG_Deseq2.txt

(cat Data/deseq_head.txt && for gene in `cat Data/salicylic_Acid_Pathway_Genes.csv | cut -f 1` ; do grep $gene FeatureCounts/*DESeq2*; done  )  > Data/SA_Deseq2.txt
 

Autophagy related genes

Heatmap Autophagy related genes transcriptionally induce/repress by 7 different viruses

DeseqOut <- read.delim("Data/ATG_Deseq2.txt", check.names = FALSE)
GeneNames <-  read.delim("Data/atg_genes_DH.txt", header = F)

colnames(GeneNames) <- c("gene","gene_name")

DeseqOut_clean <- DeseqOut %>%
  mutate(id = .[[1]]) %>%         
  select(-1) %>%
  separate(id, into = c("file", "gene"), sep = ":", extra = "merge", fill = "right") %>%
  mutate(
    file = basename(file),
    file = str_remove(file, "\\.tsv$"),
    Project = str_extract(file, "^[^_]+"),
    Condition =  gsub(".*_","",file),  
    Virus = str_match(file, "DESeq2_([^_]+)_vs_")[,2],
    sig = case_when(
      !is.na(padj) & padj < 0.001 ~ "***",
      !is.na(padj) & padj < 0.01  ~ "**",
      !is.na(padj) & padj < 0.05  ~ "*",
      TRUE ~ "ns"
    )
  )

DeseqOut_clean <- merge(GeneNames,DeseqOut_clean, by.x = "gene" )

DeseqOut_clean$Virus <- factor(DeseqOut_clean$Virus, levels = c("TuYV","CaMV","TCV","TVCV","TuMV","CMV","CaLCuV"))

Heatmap from Deseq2

#DeseqOut_clean$ColKey <- paste(DeseqOut_clean$Project,DeseqOut_clean$Condition, DeseqOut_clean$Virus, sep = "_")
DeseqOut_clean$ColKey <- DeseqOut_clean$Virus
# log2FC matrix (rows = ID (PR name), cols = Project_Infection)
log2fc_mat <- xtabs(log2FoldChange ~ gene_name + ColKey, data = DeseqOut_clean)
log2fc_mat <- as.matrix(log2fc_mat)

# helper for significance stars: take first non-empty star, else ""
first_star <- function(x) {
  x <- x[!is.na(x)]
  x <- x[x != ""]
  x <- x[!x %in% c("ns","NS","n.s.","N.S.","NA")]
  if (length(x) == 0) return("")
  x[1]
}


sig_mat <- tapply(DeseqOut_clean$sig,
                  list(DeseqOut_clean$gene_name, DeseqOut_clean$ColKey),
                  first_star)

#keep <- apply(sig_mat, 1, function(x) any(x != "" & !is.na(x)))

#sig_mat   <- sig_mat[keep, , drop = FALSE]
#log2fc_mat <- log2fc_mat[keep, , drop = FALSE]   # keep matrices aligned

# ---- colors for unscaled log2FC ----
mx <- max(abs(log2fc_mat), na.rm = TRUE)
col_fun <- colorRamp2(c(-mx, 0, mx), c("#2b8cbe", "white", "#de2d26"))
# ---- draw heatmap with stars overlaid on cells ----

colOrder <- c( "CaLCuV","TuYV","TVCV","TuMV","CMV","TCV","CaMV")
col_clusters <- c(1,2,2,1,2,2,1) 

rowOrder <- c(
  ## cluster 1
  "ATG1a","ATG7","ATG8a","ATG8b","ATG8d","ATG8g","ATG8h", "ATG13b","ATG14a","NBR1",
  ## cluster 2
  "ATG1b","ATG2","ATG3","ATG4b","ATG5","ATG6","ATG8c","ATG8e","ATG8i","ATG9","ATG10","ATG12a","ATG12b", 
  "ATG14b","ATG16", "ATG18a","ATG18f","ATG18h","ATG101",

  ## cluster 3
  "ATG1c","ATG4a","ATG8f","ATG11","ATG13a","ATG18b", "ATG18c", "ATG18d","ATG18e", "ATG18g"
)

sig_mat <- sig_mat[rowOrder,]
row_clusters <- c(rep(1,10), rep(2,19),rep(3,10)) 
ht <- Heatmap(
  log2fc_mat[rowOrder,],
  name = "log2FC",
  col = col_fun,
  na_col = "grey90",
  cluster_rows = F,
  cluster_columns = F,
  column_title = NULL, 
  row_title = NULL, 
  column_split = col_clusters, 
  column_order = colOrder, 
  row_split = row_clusters, 
  row_order = rowOrder, 
  rect_gp = gpar(col = "grey80", lwd = 0.5),
  #  stars on top of each cell
  cell_fun = function(j, i, x, y, w, h, fill) {
    s <- sig_mat[i, j]
    if (!is.na(s) && nzchar(s)) {
      grid.text(s, x, y, gp = gpar(fontsize = 12, col = "black"))
    }
  }
)

svg(filename = "Figures/CompTranscriptome_ATG_sig.svg",   width = 6, height = 8)
ht
dev.off()

pdf(file = "Figures/CompTranscriptome_ATG_sig.pdf",   width = 6, height = 8)
ht
dev.off()

SA related genes

Heatmap SA related genes transcriptionally induce/repress by 7 different viruses

DeseqOut <- read.delim("Data/SA_Deseq2.txt", check.names = FALSE)
GeneNames <-  read.delim("Data/salicylic_Acid_Pathway_Genes.csv", header = F)

colnames(GeneNames) <- c("gene","gene_name")

DeseqOut_clean <- DeseqOut %>%
  mutate(id = .[[1]]) %>%         
  select(-1) %>%
  separate(id, into = c("file", "gene"), sep = ":", extra = "merge", fill = "right") %>%
  mutate(
    file = basename(file),
    file = str_remove(file, "\\.tsv$"),
    Project = str_extract(file, "^[^_]+"),
    Condition =  gsub(".*_","",file),  
    Virus = str_match(file, "DESeq2_([^_]+)_vs_")[,2],
    sig = case_when(
      !is.na(padj) & padj < 0.001 ~ "***",
      !is.na(padj) & padj < 0.01  ~ "**",
      !is.na(padj) & padj < 0.05  ~ "*",
      TRUE ~ "ns"
    )
  )

DeseqOut_clean <- merge(GeneNames,DeseqOut_clean, by.x = "gene" )

DeseqOut_clean$Virus <- factor(DeseqOut_clean$Virus, levels = c("TuYV","CaMV","TCV","TVCV","TuMV","CMV","CaLCuV"))

Heatmap from Deseq2

DeseqOut_clean$ColKey <- DeseqOut_clean$Virus

log2fc_mat <- xtabs(log2FoldChange ~ gene_name + ColKey, data = DeseqOut_clean)
log2fc_mat <- as.matrix(log2fc_mat)


first_star <- function(x) {
  x <- x[!is.na(x)]
  x <- x[x != ""]
  x <- x[!x %in% c("ns","NS","n.s.","N.S.","NA")]
  if (length(x) == 0) return("")
  x[1]
}


sig_mat <- tapply(DeseqOut_clean$sig,
                  list(DeseqOut_clean$gene_name, DeseqOut_clean$ColKey),
                  first_star)


mx <- max(abs(log2fc_mat), na.rm = TRUE)
col_fun <- colorRamp2(c(-mx, 0, mx), c("#2b8cbe", "white", "#de2d26"))
colOrder <- c( "CaLCuV","TuYV","TVCV","TuMV","CMV","TCV","CaMV")
col_clusters <- c(1,2,2,1,2,2,1) 

rowOrder <- c(## cluster 1 
  "PR1","PR2", "FMO1",  "EDS16" ,"ICS1" ,  "PBS3",  "PR5", "ALD1" , 
   ## cluster 2
    "EDS5",  "SARD1", "CBP60g" ,  "PAD4" , "EDS1",  "NPR3", "WRKY28", 
       "PAL1"  ,  "WRKY70",   "NPR1", "NPR4", "PAL2"  ,  "MES1" ,  "MES9")

sig_mat <- sig_mat[rowOrder,]
row_clusters <- c(rep(1,8), rep(2,14)) 
ht <- Heatmap(
  log2fc_mat[rowOrder,],
  name = "log2FC",
  col = col_fun,
  na_col = "grey90",
  cluster_rows = F,
  cluster_columns = F,
  column_title = NULL, 
  row_title = NULL, 
  column_split = col_clusters, 
  column_order = colOrder, 
  row_split = row_clusters, 
 row_order = rowOrder, 
  rect_gp = gpar(col = "grey80", lwd = 0.5),
  #  stars on top of each cell
  cell_fun = function(j, i, x, y, w, h, fill) {
    s <- sig_mat[i, j]
    if (!is.na(s) && nzchar(s)) {
      grid.text(s, x, y, gp = gpar(fontsize = 12, col = "black"))
    }
  }
)

svg(filename = "Figures/CompTranscriptome_SA_sig.svg",   width = 6, height = 5)
ht
dev.off()

pdf(file = "Figures/CompTranscriptome_SA_sig.pdf",   width =5, height = 5)
ht
dev.off()

References

1. Chesnais Q, Golyaev V, Velt A, Rustenholz C, Brault V, Pooggin MM, et al. Comparative plant transcriptome profiling of arabidopsis thaliana col-0 and camelina sativa var. Celine infested with myzus persicae aphids acquiring circulative and noncirculative viruses reveals virus-and plant-specific alterations relevant to aphid feeding behavior and transmission. Microbiology spectrum. 2022;10:e00136–22.

2. Gyula P, Tóth T, Gorcsa T, Nyikó T, Sós-Hegedűs A, Szittya G. Ecotype-specific blockage of tasiARF production by two different RNA viruses in arabidopsis. Plos one. 2022;17:e0275588.

3. Wu C, Li X, Guo S, Wong S-M. Analyses of RNA-seq and sRNA-seq data reveal a complex network of anti-viral defense in TCV-infected arabidopsis thaliana. Scientific reports. 2016;6:36007.

4. Liu J-H, Lin Y, Li Y-X, Lang Z, Zhang Z, Duan C-G. A mutually antagonistic mechanism mediated by RNA m6A modification in plant-virus interactions. Nature Communications. 2025;16:10378.

5. Zorzatto C, Machado JPB, Lopes KV, Nascimento KJ, Pereira WA, Brustolini OJ, et al. NIK1-mediated translation suppression functions as a plant antiviral immunity mechanism. Nature. 2015;520:679–82.