OSCR

Molecular Drivers of Mutualistic Association Between Anemone and Anemonefish.

Code ↔ Paper

4 matches between paragraphs of the paper and lines of its authors' code, computed by the harvester (lexical-v1). Click a colored paragraph or line to see its counterpart.

The 4 matches
  1. [1] § Materials and Methods › Differential Expression Analyses ↔ clarkii_brain_DE.R, lines 389–435 · score 0.81 · clustered heatmaps, Brain_region, rlog transformed, clarkii brain, pheatmap, outlier
  2. [2] § Materials and Methods › Differential Expression Analyses ↔ Skin_DE_analysis.R, lines 524–605 · score 0.66 · log2 fold change, clarkii skin, FDR, enriched, enrichment, DE gene
  3. [3] § Materials and Methods › Sample Collection and Experimental Design ↔ clarkii_brain_DE.R, lines 389–435 · score 0.64 · optic tectum, brain stem, brain regions, cerebellum, diencephalon, telencephalon
  4. [4] § Results › Transcriptomic Responses in A. clarkii Skin Upon Acclimation ↔ Skin_DE_analysis.R, lines 524–605 · score 0.61 · log2 fold change, clarkii skin, FDR, enriched, enrichment, DE genes

Paper

Loaded from Europe PMC by your browser, not stored by OSCR: doi.org · Europe PMC

The paper is loaded when this pane is shown.

The authors' code

R · 803 lines · 32 KB · GPL-2.0 · 2 matches

  1. setwd("D:/A.clarkii/clarkii_brain_DE/featurects_output")
  3. #read featurecounts output
  4. F1_BS <- read.table("F1-BSfeaturecounts_output.txt",header=T) # Read the counts file obtained from featurecounts
  5. F1_BS_1 <- F1_BS[,c(1,7)] # Select the columns of interest - for me it was the 1st and 7th column which is gene ID and the column with the raw counts
  7. ### Repeat for all samples
  9. F2_BS <- read.table("F2-BSfeaturecounts_output.txt",header=T)
  10. F2_BS_1 <- F2_BS[,c(1,7)]
  11. F3_BS <- read.table("F3-BSfeaturecounts_output.txt",header=T)
  12. F3_BS_1 <- F3_BS[,c(1,7)]
  13. F4_BS <- read.table("F4-BSfeaturecounts_output.txt",header=T)
  14. F4_BS_1 <- F4_BS[,c(1,7)]
  15. F5_BS <- read.table("F5-BSfeaturecounts_output.txt",header=T)
  16. F5_BS_1 <- F5_BS[,c(1,7)]
  17. F6_BS <- read.table("F6-BSfeaturecounts_output.txt",header=T)
  18. F6_BS_1 <- F6_BS[,c(1,7)]
  19. F7_BS <- read.table("F7-BSfeaturecounts_output.txt",header=T)
  20. F7_BS_1 <- F7_BS[,c(1,7)]
  21. F8_BS <- read.table("F8-BSfeaturecounts_output.txt",header=T)
  22. F8_BS_1 <- F8_BS[,c(1,7)]
  23. F9_BS <- read.table("F9-BSfeaturecounts_output.txt",header=T)
  24. F9_BS_1 <- F9_BS[,c(1,7)]
  26. F1_DE <- read.table("F1-DEfeaturecounts_output.txt",header=T)
  27. F1_DE_1 <- F1_DE[,c(1,7)]
  28. F2_DE <- read.table("F2-DEfeaturecounts_output.txt",header=T)
  29. F2_DE_1 <- F2_DE[,c(1,7)]
  30. F3_DE <- read.table("F3-DEfeaturecounts_output.txt",header=T)
  31. F3_DE_1 <- F3_DE[,c(1,7)]
  32. F4_DE <- read.table("F4-DEfeaturecounts_output.txt",header=T)
  33. F4_DE_1 <- F4_DE[,c(1,7)]
  34. F5_DE <- read.table("F5-DEfeaturecounts_output.txt",header=T)
  35. F5_DE_1 <- F5_DE[,c(1,7)]
  36. F6_DE <- read.table("F6-DEfeaturecounts_output.txt",header=T)
  37. F6_DE_1 <- F6_DE[,c(1,7)]
  38. F7_DE <- read.table("F7-DEfeaturecounts_output.txt",header=T)
  39. F7_DE_1 <- F7_DE[,c(1,7)]
  40. F8_DE <- read.table("F8-DEfeaturecounts_output.txt",header=T)
  41. F8_DE_1 <- F8_DE[,c(1,7)]
  42. F9_DE <- read.table("F9-DEfeaturecounts_output.txt",header=T)
  43. F9_DE_1 <- F9_DE[,c(1,7)]
  45. F1_OT <- read.table("F1-OTfeaturecounts_output.txt",header=T)
  46. F1_OT_1 <- F1_OT[,c(1,7)]
  47. F2_OT <- read.table("F2-OTfeaturecounts_output.txt",header=T)
  48. F2_OT_1 <- F2_OT[,c(1,7)]
  49. F3_OT <- read.table("F3-OTfeaturecounts_output.txt",header=T)
  50. F3_OT_1 <- F3_OT[,c(1,7)]
  51. F4_OT <- read.table("F4-OTfeaturecounts_output.txt",header=T)
  52. F4_OT_1 <- F4_OT[,c(1,7)]
  53. F5_OT <- read.table("F5-OTfeaturecounts_output.txt",header=T)
  54. F5_OT_1 <- F5_OT[,c(1,7)]
  55. F6_OT <- read.table("F6-OTfeaturecounts_output.txt",header=T)
  56. F6_OT_1 <- F6_OT[,c(1,7)]
  57. F7_OT <- read.table("F7-OTfeaturecounts_output.txt",header=T)
  58. F7_OT_1 <- F7_OT[,c(1,7)]
  59. F8_OT <- read.table("F8-OTfeaturecounts_output.txt",header=T)
  60. F8_OT_1 <- F8_OT[,c(1,7)]
  61. F9_OT <- read.table("F9-OTfeaturecounts_output.txt",header=T)
  62. F9_OT_1 <- F9_OT[,c(1,7)]
  64. F1_TE <- read.table("F1-TEfeaturecounts_output.txt",header=T)
  65. F1_TE_1 <- F1_TE[,c(1,7)]
  66. F2_TE <- read.table("F2-TEfeaturecounts_output.txt",header=T)
  67. F2_TE_1 <- F2_TE[,c(1,7)]
  68. F3_TE <- read.table("F3-TEfeaturecounts_output.txt",header=T)
  69. F3_TE_1 <- F3_TE[,c(1,7)]
  70. F4_TE <- read.table("F4-TEfeaturecounts_output.txt",header=T)
  71. F4_TE_1 <- F4_TE[,c(1,7)]
  72. F5_TE <- read.table("F5-TEfeaturecounts_output.txt",header=T)
  73. F5_TE_1 <- F5_TE[,c(1,7)]
  74. F6_TE <- read.table("F6-TEfeaturecounts_output.txt",header=T)
  75. F6_TE_1 <- F6_TE[,c(1,7)]
  76. F7_TE <- read.table("F7-TEfeaturecounts_output.txt",header=T)
  77. F7_TE_1 <- F7_TE[,c(1,7)]
  78. F8_TE <- read.table("F8-TEfeaturecounts_output.txt",header=T)
  79. F8_TE_1 <- F8_TE[,c(1,7)]
  80. F9_TE <- read.table("F9-TEfeaturecounts_output.txt",header=T)
  81. F9_TE_1 <- F9_TE[,c(1,7)]
  83. F2_CB <- read.table("F2-CBfeaturecounts_output.txt",header=T)
  84. F2_CB_1 <- F2_CB[,c(1,7)]
  85. F3_CB <- read.table("F3-CBfeaturecounts_output.txt",header=T)
  86. F3_CB_1 <- F3_CB[,c(1,7)]
  87. F4_CB <- read.table("F4-CBfeaturecounts_output.txt",header=T)
  88. F4_CB_1 <- F4_CB[,c(1,7)]
  89. F5_CB <- read.table("F5-CBfeaturecounts_output.txt",header=T)
  90. F5_CB_1 <- F5_CB[,c(1,7)]
  91. F6_CB <- read.table("F6-CBfeaturecounts_output.txt",header=T)
  92. F6_CB_1 <- F6_CB[,c(1,7)]
  93. F7_CB <- read.table("F7-CBfeaturecounts_output.txt",header=T)
  94. F7_CB_1 <- F7_CB[,c(1,7)]
  95. F8_CB <- read.table("F8-CBfeaturecounts_output.txt",header=T)
  96. F8_CB_1 <- F8_CB[,c(1,7)]
  97. F9_CB <- read.table("F9-CBfeaturecounts_output.txt",header=T)
  98. F9_CB_1 <- F9_CB[,c(1,7)]
  100. df_list <- list(F1_BS_1,F1_DE_1,F1_OT_1,F1_TE_1,F2_BS_1,F2_DE_1,F2_OT_1,F2_TE_1,F2_CB_1,
  101. F3_BS_1,F3_DE_1,F3_OT_1,F3_TE_1,F3_CB_1,F4_BS_1,F4_DE_1,F4_OT_1,F4_TE_1,F4_CB_1,
  102. F5_BS_1,F5_DE_1,F5_OT_1,F5_TE_1,F5_CB_1,F6_BS_1,F6_DE_1,F6_OT_1,F6_TE_1,F6_CB_1,
  103. F7_BS_1,F7_DE_1,F7_OT_1,F7_TE_1,F7_CB_1,F8_BS_1,F8_DE_1,F8_OT_1,F8_TE_1,F8_CB_1,
  104. F9_BS_1,F9_DE_1,F9_OT_1,F9_TE_1,F9_CB_1)
  106. merge <- Reduce(function(x, y) merge(x, y, all=TRUE), df_list)
  107. write.csv(merge,file="clarkii_brain_cts.csv",row.names=F)
  110. #****************************************************************************
  111. #PCA from rlog transformed counts
  112. # link to guide: https://tavareshugo.github.io/data-carpentry-rnaseq/03_rnaseq_pca.html
  113. library(DESeq2)
  114. library(tidyverse)
  116. setwd("D:/A.clarkii/clarkii_brain_DE")
  117. cts <- read.csv("clarkii_brain_cts.csv",row.names="Geneid")
  118. coldata <- read.csv("brain_metadata.csv",row.names=1)
  120. #convert decimal to int
  121. cts1=round(cts)
  123. cts2 <- cts1[,-40]
  124. coldata1 <- coldata[-40,]
  127. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  128. (all(rownames(coldata1) %in% colnames(cts2)) || all(colnames(cts2) %in% rownames(coldata1)))
  130. # check if columns of count matrix are in same order as rows of sample info
  131. (all(colnames(cts2) == rownames(coldata1)))
  133. # all variables in the design formula should be converted to factors
  134. coldata1$Condition <- factor(coldata1$Condition)
  135. coldata1$Sex <- factor(coldata1$Sex)
  136. coldata1$Brain_region <- factor(coldata1$Brain_region)
  139. dds <- DESeqDataSetFromMatrix(countData = cts2, colData = coldata1, design = ~Sex + Brain_region + Condition)
  141. # Keep only genes that have non-zero reads in total
  142. keep <- rowSums(counts(dds)) > 0
  143. dds <- dds[keep,]
  145. #The rlog function transforms the count data to the log2 scale in a way which minimizes differences between samples for rows with small counts, and which normalizes with respect to library size.
  146. #The transformation is useful when checking for outliers.
  147. #Note that neither rlog transformation nor the VST are used by the differential expression estimation in DESeq, which always occurs on the raw count data, through generalized linear modeling which incorporates knowledge of the variance-mean dependence. The rlog transformation and VST are offered as separate functionality which can be used for visualization, clustering or other machine learning tasks.
  149. rlogcounts <- assay(rlog(dds, blind=TRUE)) #blind = T will make the log transformation blind to the sample info)
  151. #Use the rlog transformed counts to do a PCA
  152. #To use prcomp, the samples should be the rows and the variables (in this case genes) should be the columns. So we need to transpose the data
  153. transpose <- t(rlogcounts)
  155. sample_pca <- prcomp(transpose,center=T, scale=T)#Often it is a good idea to standardize the variables before doing the PCA. This is often done by centering the data on the mean and then scaling it by dividing by the standard deviation. This ensures that the PCA is not too influenced by genes with higher absolute expression. By default, the prcomp() function does the centering but not the scaling.
  157. pc_eigenvalues <- sample_pca$sdev^2 # to extract the variance explained by each PC from our sample_pca object.
  159. #pc_eigenvalues is a vector. So to make a plot with ggplot2, we need to transform it to a dataframe/tibble
  161. # create a "tibble" manually with
  162. # a variable indicating the PC number
  163. # and a variable with the variances
  164. pc_eigenvalues <- tibble(PC = factor(1:length(pc_eigenvalues)),
  165. variance = pc_eigenvalues) %>%
  166. # add a new column with the percent variance
  167. mutate(pct = variance/sum(variance)*100) %>%
  168. # add another column with the cumulative variance explained
  169. mutate(pct_cum = cumsum(pct))
  171. # print the result
  172. pc_eigenvalues
  174. # The pc_eigenvalues table can now be used to produce a Scree Plot, which shows the fraction of total variance explained by each principal component.
  175. # The plot will show both the variance explained by individual PCs as well as the cumulative variance, using a type of visualisation known as a pareto chart
  176. #A Pareto chart is a type of chart that contains both bars and a line graph, where individual values are represented in descending order by bars, and the cumulative total is represented by the line.
  178. p <- pc_eigenvalues %>%
  179. ggplot(aes(x = PC)) +
  180. geom_col(aes(y = pct)) +
  181. geom_line(aes(y = pct_cum, group = 1)) +
  182. geom_point(aes(y = pct_cum)) +
  183. labs(x = "Principal component", y = "Fraction variance explained")
  184. p
  185. # The plot shows how successive PCs explain less and less of the total variance in the original data. Also note that in this case 7 components are enough to virtually explain all of the variance in the dataset. This makes sense in this case since we only have 7 biological samples.
  187. p <- p + theme_bw() + theme(axis.text=element_text(size=26), axis.title=element_text(size=30,face = "bold"))
  188. p
  190. setwd("./figures")
  191. ggsave("PCA_variance_for_each_PC_remvedF9BS.png", dpi = 300,units="cm",width=60,height=40)
  193. #Visualising samples on PC space
  194. library(ggpubr)
  195. setwd("D:/A.clarkii/clarkii_brain_DE")
  196. #extract the PC score for each sample from the prcomp object
  197. coldata <- read.csv("brain_metadata.csv")
  198. coldata1 <- coldata[-40,]
  200. pc_scores <- sample_pca$x ## The PC scores are stored in the "x" value of the prcomp object
  202. #The pc_scores object is of class matrix, so we need to first convert it to a data.frame/tibble for ggplot2.
  203. pc_scores <- pc_scores %>%
  204. # convert to a tibble retaining the sample names as a new column
  205. as_tibble(rownames = "Sample")
  207. # print the result
  208. pc_scores
  210. library(ggrepel)
  211. # Finally plot the PCA
  212. p <- pc_scores %>%
  213. # join with "sample_info" table
  214. full_join(coldata1, by = "Sample") %>%
  215. # create the plot
  216. ggplot(aes(x = PC1, y = PC2, colour = Brain_region,label=Sample,group=Brain_region)) +
  217. geom_point(size=5,aes(shape=interaction(Condition,Sex))) + scale_shape_manual(values = c(1, 16, 2, 17)) +
  218. geom_text_repel(size=5,max.overlaps = 15,min.segment.length = 0) +
  219. stat_conf_ellipse(bary=F, level = 0.95,npoint = 1000, geom = "polygon",alpha=0.3, aes (fill=Brain_region))
  221. #library(scales) - to see what default colors are used
  222. #show_col(hue_pal()(5))
  224. p
  226. #no labels
  227. p <- pc_scores %>%
  228. # join with "sample_info" table
  229. full_join(coldata, by = "Sample") %>%
  230. # create the plot
  231. ggplot(aes(x = PC1, y = PC2, colour = Brain_region,group=Brain_region)) +
  232. geom_point(size=5,aes(shape=interaction(Condition,Sex))) + scale_shape_manual(values = c(1, 16, 2, 17)) +
  233. stat_conf_ellipse(bary=F, level = 0.95,npoint = 1000, geom = "polygon",alpha=0.3, aes (fill=Brain_region))
  236. pc1_variance <- pc_eigenvalues[1,3]
  237. pc2_variance <- pc_eigenvalues[2,3]
  239. p <- p + xlab(paste0("PC1: ",round(pc1_variance,digits = 2),"%")) +
  240. ylab(paste0("PC2: ",round(pc2_variance,digits = 2),"%"))
  242. p
  244. p <- p + theme_bw() + theme(axis.text=element_text(size=26), axis.title=element_text(size=30,face = "bold"))
  245. p
  247. p <- p + theme(legend.title = element_text(size=30),
  248. legend.text = element_text(size=26))
  250. p
  252. setwd("./figures")
  253. ggsave("PCA_brain_all_factors_noLabels_removedF9BS.png", dpi = 300)
  255. #Exploring correlation between genes and PCs
  256. #This is to see which genes have the most influence on each PC axis. This information is contained in the variable loadings of the PCA, within the rotation value of the prcomp object.
  258. pc_loadings <- sample_pca$rotation
  260. #The pc_loading object is a matrix, which we can convert to a data.frame/tibble for plotting and data manipulation
  262. pc_loadings <- pc_loadings %>%
  263. as_tibble(rownames = "gene")
  265. # print the result
  266. pc_loadings
  268. #extract the top 10 genes with highest loading on PC1 and PC2
  269. top_genes <- pc_loadings %>%
  270. # select only the PCs we are interested in
  271. select(gene, PC1, PC2) %>%
  272. # convert to a "long" format
  273. pivot_longer(matches("PC"), names_to = "PC", values_to = "loading") %>%
  274. # for each PC
  275. group_by(PC) %>%
  276. # arrange by descending order of loading
  277. arrange(desc(abs(loading))) %>%
  278. # take the 10 top rows
  279. slice(1:10) %>%
  280. # pull the gene column as a vector
  281. pull(gene) %>%
  282. # ensure only unique genes are retained
  283. unique()
  285. top_genes
  287. #Now, we can use this list of gene names to subset the eigenvalues table:
  288. top_loadings <- pc_loadings %>%
  289. filter(gene %in% top_genes)
  291. #plot
  292. loadings_plot <- ggplot(data = top_loadings) +
  293. geom_segment(aes(x = 0, y = 0, xend = PC1, yend = PC2),
  294. arrow = arrow(length = unit(0.1, "in")),
  295. colour = "black",linewidth=1) +
  296. geom_text_repel(aes(x = PC1, y = PC2, label = gene),
  297. size = 8,hjust=2, vjust=3) + theme_bw() +
  298. theme(axis.text=element_text(size=26),
  299. axis.title=element_text(size=30,face = "bold")) +
  300. xlab("PC1") + ylab ("PC2")
  302. loadings_plot
  304. #Interpretation
  305. #Use the loading plot to identify which variables have the largest effect on each component.
  306. #Loadings can range from -1 to 1. Loadings close to -1 or 1 indicate that the variable strongly influences the component.
  307. #Loadings close to 0 indicate that the variable has a weak influence on the component.
  308. #Positive loadings indicate a variable and a principal component are positively correlated. Negative loadings indicate a negative correlation.
  310. setwd("./figures")
  311. ggsave("PCA_loadings_top10genesForeachPC_RemovedF9BS.png", dpi = 300)
  313. # Combining PC scores and PC loadings to see correlation - biplot
  314. #For example, geneA, with a negative score of -0.00605 (loading score), will be associated with whatever samples are shifted toward the negative end of your PC1 on a bi-plot.
  316. if (!requireNamespace('BiocManager', quietly = TRUE))
  317. install.packages('BiocManager')
  319. BiocManager::install('PCAtools')
  321. library(PCAtools)
  322. setwd("D:/A.clarkii/clarkii_brain_DE")
  324. cts <- read.csv("clarkii_brain_cts.csv",row.names="Geneid")
  325. coldata <- read.csv("brain_metadata.csv",row.names = 1)
  327. #convert decimal to int
  328. cts1=round(cts)
  330. cts2 <- cts1[,-40]
  331. coldata1 <- coldata[-40,]
  334. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  335. (all(rownames(coldata1) %in% colnames(cts2)) || all(colnames(cts2) %in% rownames(coldata1)))
  337. # check if columns of count matrix are in same order as rows of sample info
  338. (all(colnames(cts2) == rownames(coldata1)))
  340. # all variables in the design formula should be converted to factors
  341. coldata1$Condition <- factor(coldata1$Condition)
  342. coldata1$Sex <- factor(coldata1$Sex)
  343. coldata1$Brain_region <- factor(coldata1$Brain_region)
  347. dds <- DESeqDataSetFromMatrix(countData = cts2, colData = coldata1, design = ~Sex + Condition + Brain_region)
  349. # Keep only genes that have non-zero reads in total
  350. keep <- rowSums(counts(dds)) > 0
  351. dds <- dds[keep,]
  353. #The rlog function transforms the count data to the log2 scale in a way which minimizes differences between samples for rows with small counts, and which normalizes with respect to library size.
  354. #The transformation is useful when checking for outliers.
  355. #Note that neither rlog transformation nor the VST are used by the differential expression estimation in DESeq, which always occurs on the raw count data, through generalized linear modeling which incorporates knowledge of the variance-mean dependence. The rlog transformation and VST are offered as separate functionality which can be used for visualization, clustering or other machine learning tasks.
  357. rlogcounts <- assay(rlog(dds, blind=TRUE)) #blind = T will make the log transformation blind to the sample info)
  360. a <- pca(rlogcounts, metadata = coldata1,center=T,scale=T)
  362. biplot(a, showLoadings = TRUE, colby='Brain_region', shape='Sex',
  363. lab = NULL, pointSize = 5, sizeLoadingsNames = 7,ellipse=F,
  364. ntopLoadings = 10,max.overlaps = 30, boxedLoadingsNames = FALSE,
  365. legendPosition="right",legendLabSize = 20,axisLabSize=26,
  366. legendIconSize = 8,colLegendTitle=NULL,shapeLegendTitle=NULL)
  368. setwd("./figures")
  369. ggsave("PCA_biplot_removedF9BS.png", dpi = 300)
  372. #****************************************************************************
  373. #clustering heatmap to check if samples group by sex
  374. library(DESeq2)
  375. library("pheatmap")
  376. library("dplyr")
  378. setwd("D:/A.clarkii/clarkii_brain_DE")
  380. cts <- read.csv("clarkii_brain_cts.csv",row.names="Geneid")
  381. coldata <- read.csv("brain_metadata.csv",row.names = 1)
  383. #convert decimal to int
  384. cts1=round(cts)
  386. cts2 <- data.frame(cts1[,c(-40)])
  387. coldata1 <- coldata[c(-40),]
  389. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  390. (all(rownames(coldata1) %in% colnames(cts2)) || all(colnames(cts2) %in% rownames(coldata1)))
  392. # check if columns of count matrix are in same order as rows of sample info
  393. (all(colnames(cts2) == rownames(coldata1)))
  396. # all variables in the design formula should be converted to factors
  397. coldata1$Condition <- factor(coldata1$Condition)
  398. coldata1$Sex <- factor(coldata1$Sex)
  399. coldata1$Brain_region <- factor(coldata1$Brain_region)
  402. dds <- DESeqDataSetFromMatrix(countData = cts2, colData = coldata1, design = ~Condition + Sex + Brain_region)
  404. # Keep only genes that have non-zero reads in total
  405. keep <- rowSums(counts(dds)) > 0
  406. dds <- dds[keep,]
  408. #The rlog function transforms the count data to the log2 scale in a way which minimizes differences between samples for rows with small counts, and which normalizes with respect to library size.
  409. #The transformation is useful when checking for outliers.
  410. #Note that neither rlog transformation nor the VST are used by the differential expression estimation in DESeq, which always occurs on the raw count data, through generalized linear modeling which incorporates knowledge of the variance-mean dependence. The rlog transformation and VST are offered as separate functionality which can be used for visualization, clustering or other machine learning tasks.
  412. rlogcounts <- assay(rlog(dds, blind=TRUE)) #blind = T will make the log transformation blind to the sample info)
  414. cor <- cor(rlogcounts) #Calculate the correlation values between samples
  416. library(RColorBrewer)
  417. display.brewer.all(colorblindFriendly = T)
  418. heat_colors <- brewer.pal(9, "Greys")
  420. annotation_cols <- list(Condition = c(Control = "#808080", Interaction = "#FEBE00"),
  421. Sex = c(Male = "#6B8E23", Female = "#CC8899"),
  422. Brain_region = c("Brain stem" = "#5BC8E8",
  423. "Cerebellum" = "#1FC4B0",
  424. "Diencephalon" = "#8E9CE0",
  425. "Optic tectum" = "#A0A018",
  426. "Telencephalon" = "#F26FD0"))
  429. p <- pheatmap(cor,color=heat_colors, cellheight=35,cellwidth=40,cluster_rows=TRUE, show_rownames=T,show_colnames=T,legend=T,fontsize = 30,
  430. cluster_cols=T, annotation=select(coldata1,c(Condition,Sex,Brain_region)),,annotation_colors = annotation_cols)
  432. setwd("./figures/")
  433. png("clustering_heatmap_removedF9BS.png",res = 300,units="cm", width = 82, height = 63)
  434. print(p)
  435. dev.off()
  438. #****************************************************************************
  439. #LRT to see the effect of sex, brain_region and condition on gene expression
  440. library(DESeq2)
  442. setwd("D:/A.clarkii/clarkii_brain_DE")
  444. cts <- read.csv("clarkii_brain_cts.csv",row.names="Geneid")
  445. coldata <- read.csv("brain_metadata.csv",row.names = 1)
  447. #convert decimal to int
  448. cts1=round(cts)
  450. cts2 <- data.frame(cts1[,c(-40)])
  451. coldata1 <- coldata[c(-40),]
  453. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  454. (all(rownames(coldata1) %in% colnames(cts2)) || all(colnames(cts2) %in% rownames(coldata1)))
  456. # check if columns of count matrix are in same order as rows of sample info
  457. (all(colnames(cts2) == rownames(coldata1)))
  460. # all variables in the design formula should be converted to factors
  461. coldata1$Condition <- factor(coldata1$Condition)
  462. coldata1$Sex <- factor(coldata1$Sex)
  463. coldata1$Brain_region <- factor(coldata1$Brain_region)
  466. dds <- DESeqDataSetFromMatrix(countData = cts2, colData = coldata1,
  467. design = ~Sex + Condition + Brain_region)
  469. # Keep only genes that have non-zero reads in total
  470. keep <- rowSums(counts(dds)) > 0
  471. dds <- dds[keep,]
  473. # total expressed genes in brain (passed the rowSums > 0 filter)
  474. n_expressed_brain <- nrow(dds)
  475. n_expressed_brain #24907
  477. # testing for sex effect
  478. dds_LRT_sex <- DESeq(dds, test = "LRT", reduced = ~ Condition + Brain_region)
  479. res_sex <- results(dds_LRT_sex)
  480. res_sex
  481. summary(res_sex)
  482. results <- data.frame(res_sex)
  483. results <- na.omit(results)
  485. res_sig_padj <- results[results$padj < 0.05, ]
  487. # number significant
  488. n_sig_sex <- sum(res_sex$padj < 0.05, na.rm = TRUE)
  489. n_sig_sex #365
  491. # proportion of expressed genes
  492. round(100 * n_sig_sex / n_expressed_brain, 1)
  494. setwd("./LRT_results/")
  495. write.csv(results,file="LRT_results_sex_effect.csv")
  496. write.csv(res_sig_padj,file="LRT_results_sex_effect_sig.csv")
  498. # testing for condition effect
  499. dds_LRT_condition <- DESeq(dds, test = "LRT", reduced = ~ Sex + Brain_region)
  500. res_condition <- results(dds_LRT_condition)
  501. res_condition
  502. summary(res_condition)
  503. results <- data.frame(res_condition)
  504. results <- na.omit(results)
  506. res_sig_padj <- results[results$padj<0.05,]
  508. setwd("./LRT_results/")
  509. write.csv(results,file="LRT_results_condition_effect.csv")
  510. write.csv(res_sig_padj,file="LRT_results_condition_effect_sig.csv")
  512. # testing for Brain_region effect
  513. dds_LRT_BrainRegion <- DESeq(dds, test = "LRT", reduced = ~ Sex + Condition)
  514. res_BrainRegion <- results(dds_LRT_BrainRegion)
  515. res_BrainRegion
  516. summary(res_BrainRegion)
  517. results <- data.frame(res_BrainRegion)
  518. results <- na.omit(results)
  520. res_sig_padj <- results[results$padj<0.05,]
  522. setwd("./LRT_results/")
  523. write.csv(results,file="LRT_results_BrainRegion_effect.csv")
  524. write.csv(res_sig_padj,file="LRT_results_BrainRegion_effect_sig.csv")
  526. #****************************************************************************
  527. # DE analysis - wald test
  528. library(DESeq2)
  530. setwd("D:/A.clarkii/clarkii_brain_DE")
  532. cts <- read.csv("clarkii_brain_cts.csv",row.names="Geneid")
  533. coldata <- read.csv("brain_metadata.csv",row.names = 1)
  535. #convert decimal to int
  536. cts1=round(cts)
  538. cts2 <- data.frame(cts1[,c(-40)])
  539. coldata1 <- coldata[c(-40),]
  541. #### BS
  542. cts_BS <- cts2[,c(1,5,10,15,20,25,30,35)]
  543. coldata_BS <- coldata1[c(1,5,10,15,20,25,30,35),]
  545. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  546. (all(rownames(coldata_BS) %in% colnames(cts_BS)) || all(colnames(cts_BS) %in% rownames(coldata_BS)))
  548. # check if columns of count matrix are in same order as rows of sample info
  549. (all(colnames(cts_BS) == rownames(coldata_BS)))
  551. dds_BS <- DESeqDataSetFromMatrix(countData = cts_BS, colData = coldata_BS,design = ~ Sex + Condition)
  552. levels(dds_BS$Condition)
  553. levels(dds_BS$Sex)
  555. dds_BS$condition <- relevel(dds_BS$Condition, ref = "Control")
  558. # Keep only genes that have non-zero reads in total
  559. keep <- rowSums(counts(dds_BS)) > 0
  560. dds_BS <- dds_BS[keep,]
  562. dds_BS <- DESeq(dds_BS)
  563. resultsNames(dds_BS)
  565. resBS <- results(dds_BS, alpha = 0.05, contrast = c("Condition","Interaction","Control"), test ="Wald")#So up/downregulation refers to interaction i.e upregulated in interaction or downregulated in interaction
  566. summary(resBS)
  567. resBS
  569. resultsBS <- data.frame(resBS)
  570. resultsBS <- na.omit(resultsBS)
  571. sig_results_BS <- resultsBS[resultsBS$padj<0.05, ]
  573. setwd("./DE_results_wald_test/")
  574. write.csv(resultsBS,file="DE_results_all_BS.csv")
  575. write.csv(sig_results_BS,file="DE_results_sig_BS.csv")
  577. #### OT
  578. cts_OT <- cts2[,c(3,7,12,17,22,27,32,37,41)]
  579. coldata_OT <- coldata1[c(3,7,12,17,22,27,32,37,41),]
  581. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  582. (all(rownames(coldata_OT) %in% colnames(cts_OT)) || all(colnames(cts_OT) %in% rownames(coldata_OT)))
  584. # check if columns of count matrix are in same order as rows of sample info
  585. (all(colnames(cts_OT) == rownames(coldata_OT)))
  587. dds_OT <- DESeqDataSetFromMatrix(countData = cts_OT, colData = coldata_OT,design = ~ Sex + Condition)
  588. levels(dds_OT$Condition)
  589. levels(dds_OT$Sex)
  591. dds_OT$condition <- relevel(dds_OT$Condition, ref = "Control")
  594. # Keep only genes that have non-zero reads in total
  595. keep <- rowSums(counts(dds_OT)) > 0
  596. dds_OT <- dds_OT[keep,]
  598. dds_OT <- DESeq(dds_OT)
  599. resultsNames(dds_OT)
  601. resOT <- results(dds_OT, alpha = 0.05, contrast = c("Condition","Interaction","Control"), test ="Wald")#So up/downregulation refers to interaction i.e upregulated in interaction or downregulated in interaction
  602. summary(resOT)
  603. resOT
  605. resultsOT <- data.frame(resOT)
  606. resultsOT <- na.omit(resultsOT)
  607. sig_results_OT <- resultsOT[resultsOT$padj<0.05, ]
  609. setwd("./DE_results_wald_test/")
  610. write.csv(resultsOT,file="DE_results_all_OT.csv")
  611. write.csv(sig_results_OT,file="DE_results_sig_OT.csv")
  614. #### TE
  615. cts_TE <- cts2[,c(4,8,13,18,23,28,33,38,42)]
  616. coldata_TE <- coldata1[c(4,8,13,18,23,28,33,38,42),]
  618. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  619. (all(rownames(coldata_TE) %in% colnames(cts_TE)) || all(colnames(cts_TE) %in% rownames(coldata_TE)))
  621. # check if columns of count matrix are in same order as rows of sample info
  622. (all(colnames(cts_TE) == rownames(coldata_TE)))
  624. dds_TE <- DESeqDataSetFromMatrix(countData = cts_TE, colData = coldata_TE,design = ~ Sex + Condition)
  625. levels(dds_TE$Condition)
  626. levels(dds_TE$Sex)
  628. dds_TE$condition <- relevel(dds_TE$Condition, ref = "Control")
  631. # Keep only genes that have non-zero reads in total
  632. keep <- rowSums(counts(dds_TE)) > 0
  633. dds_TE <- dds_TE[keep,]
  635. dds_TE <- DESeq(dds_TE)
  636. resultsNames(dds_TE)
  638. resTE <- results(dds_TE, alpha = 0.05, contrast = c("Condition","Interaction","Control"), test ="Wald")#So up/downregulation refers to interaction i.e upregulated in interaction or downregulated in interaction
  639. summary(resTE)
  640. resTE
  642. resultsTE <- data.frame(resTE)
  643. resultsTE <- na.omit(resultsTE)
  644. sig_results_TE <- resultsTE[resultsTE$padj<0.05, ]
  646. setwd("./DE_results_wald_test/")
  647. write.csv(resultsTE,file="DE_results_all_TE.csv")
  648. write.csv(sig_results_TE,file="DE_results_sig_TE.csv")
  651. #### CB
  652. cts_CB <- cts2[,c(9,14,19,24,29,34,39,43)]
  653. coldata_CB <- coldata1[c(9,14,19,24,29,34,39,43),]
  655. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  656. (all(rownames(coldata_CB) %in% colnames(cts_CB)) || all(colnames(cts_CB) %in% rownames(coldata_CB)))
  658. # check if columns of count matrix are in same order as rows of sample info
  659. (all(colnames(cts_CB) == rownames(coldata_CB)))
  661. dds_CB <- DESeqDataSetFromMatrix(countData = cts_CB, colData = coldata_CB,design = ~ Sex + Condition)
  662. levels(dds_CB$Condition)
  663. levels(dds_CB$Sex)
  665. dds_CB$condition <- relevel(dds_CB$Condition, ref = "Control")
  668. # Keep only genes that have non-zero reads in total
  669. keep <- rowSums(counts(dds_CB)) > 0
  670. dds_CB <- dds_CB[keep,]
  672. dds_CB <- DESeq(dds_CB)
  673. resultsNames(dds_CB)
  675. resCB <- results(dds_CB, alpha = 0.05, contrast = c("Condition","Interaction","Control"), test ="Wald")#So up/downregulation refers to interaction i.e upregulated in interaction or downregulated in interaction
  676. summary(resCB)
  677. resCB
  679. resultsCB <- data.frame(resCB)
  680. resultsCB <- na.omit(resultsCB)
  681. sig_results_CB <- resultsCB[resultsCB$padj<0.05, ]
  683. setwd("./DE_results_wald_test/")
  684. write.csv(resultsCB,file="DE_results_all_CB.csv")
  685. write.csv(sig_results_CB,file="DE_results_sig_CB.csv")
  688. #### DE
  689. cts_DE <- cts2[,c(2,6,11,16,21,26,31,36,40)]
  690. coldata_DE <- coldata1[c(2,6,11,16,21,26,31,36,40),]
  692. # check if all colnames of count matrix are in the rownames of sample info; if false exit the script
  693. (all(rownames(coldata_DE) %in% colnames(cts_DE)) || all(colnames(cts_DE) %in% rownames(coldata_DE)))
  695. # check if columns of count matrix are in same order as rows of sample info
  696. (all(colnames(cts_DE) == rownames(coldata_DE)))
  698. dds_DE <- DESeqDataSetFromMatrix(countData = cts_DE, colData = coldata_DE,design = ~ Sex + Condition)
  699. levels(dds_DE$Condition)
  700. levels(dds_DE$Sex)
  702. dds_DE$condition <- relevel(dds_DE$Condition, ref = "Control")
  705. # Keep only genes that have non-zero reads in total
  706. keep <- rowSums(counts(dds_DE)) > 0
  707. dds_DE <- dds_DE[keep,]
  709. dds_DE <- DESeq(dds_DE)
  710. resultsNames(dds_DE)
  712. resDE <- results(dds_DE, alpha = 0.05, contrast = c("Condition","Interaction","Control"), test ="Wald")#So up/downregulation refers to interaction i.e upregulated in interaction or downregulated in interaction
  713. summary(resDE)
  714. resDE
  716. resultsDE <- data.frame(resDE)
  717. resultsDE <- na.omit(resultsDE)
  718. sig_results_DE <- resultsDE[resultsDE$padj<0.05, ]
  720. setwd("./DE_results_wald_test/")
  721. write.csv(resultsDE,file="DE_results_all_DE.csv")
  722. write.csv(sig_results_DE,file="DE_results_sig_DE.csv")
  725. #****************************************************************************
  726. # merge DE gene list with GO functional annotation
  727. setwd("D:/A.clarkii/clarkii_brain_DE/DE_results_wald_test/")
  729. GO <- read.csv("A.clarkii_GO_Annotation_edited.csv")
  731. BS <- read.csv("DE_results_sig_final_BS.csv")
  732. merge <- merge(BS,GO,by="Geneid",all.x=T)
  733. write.csv(merge,file="DE_results_sig_final_BS_annotation.csv")
  735. CB <- read.csv("DE_results_sig_final_CB.csv")
  736. merge <- merge(CB,GO,by="Geneid",all.x=T)
  737. write.csv(merge,file="DE_results_sig_final_CB_annotation.csv")
  739. DE <- read.csv("DE_results_sig_final_DE.csv")
  740. merge <- merge(DE,GO,by="Geneid",all.x=T)
  741. write.csv(merge,file="DE_results_sig_final_DE_annotation.csv")
  743. TE <- read.csv("DE_results_sig_final_TE.csv")
  744. merge <- merge(TE,GO,by="Geneid",all.x=T)
  745. write.csv(merge,file="DE_results_sig_final_TE_annotation.csv")
  747. OT <- read.csv("DE_results_sig_final_OT.csv")
  748. merge <- merge(OT,GO,by="Geneid",all.x=T)
  749. write.csv(merge,file="DE_results_sig_final_OT_annotation.csv")
  751. ############################################## Venn Diagrams
  752. # http://www.interactivenn.net/index2.html
  756. ##### merge individual brain regions to get a combined mega file
  758. setwd("D:/A.clarkii/clarkii_brain_DE/DE_results_wald_test/final_files")
  760. BS <- read.csv("DE_results_sig_final_BS_annotation.csv")
  761. CB <- read.csv("DE_results_sig_final_CB_annotation.csv")
  762. DE <- read.csv("DE_results_sig_final_DE_annotation.csv")
  763. OT <- read.csv("DE_results_sig_final_OT_annotation.csv")
  764. TE <- read.csv("DE_results_sig_final_TE_annotation.csv")
  766. df_list <- list(BS,CB,DE,OT,TE)
  768. merge <- Reduce(function(x, y) merge(x, y, all=TRUE), df_list)
  769. write.csv(merge,file="clarkii_DE_final_Allcombined.csv",row.names=F)
  771. ############################################## Heatmap
  772. library(ggplot2)
  773. library(ggbeeswarm) # for beeswarm plot
  774. library(ggrepel)
  776. setwd("D:/A.clarkii/clarkii_brain_DE/DE_results_wald_test/final_files")
  778. file <- read.csv("file_for_box.csv")
  779. file1 <- na.omit(file)
  781. ggplot(file1,aes(x=Brain_Region,y=log2FoldChange,color=Brain_Region)) +
  782. geom_quasirandom(size=4) + theme_bw() + theme(legend.title=element_blank()) +
  783. theme(legend.text=element_text(size=18)) + geom_text_repel(data=file1,aes(label=Gene_symbol_NCBI),position=position_quasirandom(),size=4) +
  784. theme(axis.text=element_text(size=18),axis.title.x = element_blank(), axis.title.y=element_text(size=22)) +
  785. scale_fill_manual(values=c("#F8766D","A3A500","00BF7D","00B0F6","E76BF3")) +
  786. geom_hline(yintercept=0, linetype="dashed",color="black",size=1)
  788. setwd("D:/A.clarkii/clarkii_brain_DE/figures")
  789. ggsave(file="brain_regions_beeswarm.png",dpi=300,width=20, height=20,units="cm")
  791. # ============================================
  792. # Test overlap between sig sex effect genes from LRT and the sig DE genes from Wald test
  793. # ============================================
  795. setwd("D:/A.clarkii/clarkii_brain_DE")
  797. lrt_sex <- read.csv("LRT_results/LRT_results_sex_effect_sig.csv",row.names=1)
  798. sex <- rownames(lrt_sex)
  800. de_genes <- read.csv("DE_results_wald_test/final_files/clarkii_DE_final_Allcombined.csv",row.names=1)
  801. brain_de <- rownames(de_genes)
  803. overlap <- intersect(sex, brain_de)

clarkii_brain_DE.R at commit da1ba09, under GPL-2.0 · at the source

Overview

Authors: Sneha Suresh1,2, Daniele Romeo1, Celia Schunter1,3
ORCID iDs: Sneha Suresh, Daniele Romeo, Celia Schunter
  1. Swire Institute of Marine Science and School of Biological Sciences, The University of Hong Kong, Hong Kong SAR, China
  2. Department of Environmental Conservation, University of Massachusetts, Amherst, MA 01003, USA
  3. State Key Laboratory of Marine Environmental Health (SKLMEH), City University of Hong Kong, Hong Kong SAR, China
Institutions: University of Massachusetts Amherst (United States); University of Hong Kong (Hong Kong SAR China); City University of Hong Kong (Hong Kong SAR China)
Journal: Genome biology and evolution, volume 18, issue 8, article evag197
Dates: accepted 24 July 2026; published online 31 July 2026; in print August 2026
Type: Research article · Language: English
License: CC BY-NC
Identifiers: DOI 10.1093/gbe/evag197 · PMID 42535633 · PMCID PMC13481981 · OpenAlex W4413614195
Open access: gold, a free copy (OpenAlex)
Preprint: doi:10.1101/2025.08.21.671668
Status: code verified
Categories: genetics / omics (modality), other (organism), cellular / molecular (subfield)
Methods: Statistics, Smoothing, state filtering, decompositions, Connectivity
Keywords: mutualism, transcriptomics, symbiosis, anemonefish, sea anemone, Amphiprion clarkii
MeSH: Anemone*, Fishes*, Gene Expression Profiling*, Symbiosis*, Animals, Aromatase, Behavior, Animal, Brain, Cnidarian Venoms, Cytoskeletal Proteins, Fish Proteins, Olfactory Perception, Skin (* major topic)
Topic: Marine Invertebrate Physiology and Ecology (Paleontology, Earth and Planetary Sciences), according to OpenAlex
Funding: University of Hong Kong
Citations: not cited yet (Europe PMC); 96 references in the paper

Abstract

The anemone–anemonefish mutualism is one of the most iconic in the marine environment. While the evolution of this mutualistic relationship has contributed to the ecological success of both partners, the underlying molecular processes that establish and maintain it remain poorly understood, particularly how anemonefish tolerate anemone venom. Here, we characterize the transcriptional dynamics in both the anemonefish Amphiprion clarkii and its host anemone Entacmaea quadricolor 48 h after association, providing a rare insight into the coordinated molecular processes in both partners that underlie symbiosis establishment. Upon acclimation with an anemone, anemonefish showed differential regulation of sensory perception and memory genes in key brain regions, indicating activation of neural pathways that may facilitate host recognition and mutualism establishment. In the fish's skin, altered expression of genes involved in neurotransmitter release, cytoskeleton organization, and venom receptor proteins points to mechanisms of resistance to anemone venom. This resistance is particularly remarkable since anemone hosting fish exhibited increased expression of genes encoding mechanoreceptors, putative venom-associated proteins, and ion channels involved in nematocyst discharge, indicating the anemone does indeed mount an active response to their mutualistic partner. By simultaneously capturing the molecular responses of both symbiotic partners, our results reveal the complex, coordinated interplay of molecular events in both species that play a pivotal role in establishing this mutualistic relationship.

Reproduced under the paper's license (CC BY-NC), from the paper cited above.

Repository

Its files are read in the Code ↔ Paper reader above, with 4 matches between paragraphs and lines of code.

sneha100895/Anemone-anemonefish-mutualism-RNA-Seq

License: GPL-2.0
State: the link answers, verified on 26 September 2026
Evidence: files inventoried
Commit: da1ba091cf1a7213c50910e8b8c86c83ba1aaadc, 28 July 2026
Languages: R (3)
Size: 6 files, 3 scripts
Software Heritage: not archived
Found in: “Data Availability”
Holds: README, license file
Not found: CITATION.cff, environment file, tests, continuous integration, documentation
Tools: DESeq2 (3 files), ggplot2 (3 files), ggpubr (3 files), pheatmap (3 files), tidyverse (3 files)
Availability: 1 check, the latest on 26 September 2026: the link answers
  • 26 September 2026: the link answers
5 files
At the source: github.com/sneha100895/Anemone-anemonefish-mutualism-RNA-Seq

The paper's code and data availability statement is in the Data section.

Tracing map

Proposed by the machine: these links were found in the paper and verified at the source, without human review. The map will receive a Zenodo DOI once one of the paper's authors has validated it with their ORCID.

What the map holds:

  • 1 repository of the authors' code, each at its verified commit, with its license and how the link was found in the paper;
  • 3 scripts, each with its path and the digest of its content;
  • 4 matches between paragraphs of the paper and lines of the code (method lexical-v1);
  • neither the text of the paper nor the code itself.

Its JSON (tracing-map.json) is deposited on Zenodo with its DOI once the map is validated.

Data

Datasets cited

Data Availability

The RNA-Seq raw sequences are deposited in NCBI under BioProject ID PRJNA1056457. Code used for analyses is available at https://github.com/sneha100895/Anemone-anemonefish-mutualism-RNA-Seq

Reproduced under the paper's license (CC BY-NC), from the paper cited above.

Versions

The history of this record: each version stored by the harvester or made by a correction of its authors or of the maintainers of its code, and what changed in its facts. The texts of the paper (its abstract, its availability statements) are not part of it; versions that changed only those are not listed.

Version 1, 27 September 2026: the first record

Recorded: type, language, journal, volume, issue, pages, dates, 3 authors, 6 keywords, 13 MeSH terms, 1 funder, 92 references.

Cite

This paper

Suresh, S., Romeo, D., & Schunter, C. (2026). Molecular Drivers of Mutualistic Association Between Anemone and Anemonefish. Genome biology and evolution, 18(8), evag197. https://doi.org/10.1093/gbe/evag197

BibTeX

@article{suresh2026molecular,
author = {Suresh, Sneha and Romeo, Daniele and Schunter, Celia},
title = {{Molecular Drivers of Mutualistic Association Between Anemone and Anemonefish}},
journal = {Genome biology and evolution},
year = {2026},
month = aug,
volume = {18},
number = {8},
pages = {evag197},
publisher = {Oxford University Press},
issn = {1759-6653},
doi = {10.1093/gbe/evag197},
url = {https://doi.org/10.1093/gbe/evag197},
pmid = {42535633},
pmcid = {PMC13481981}
}

RIS

TY - JOUR
AU - Suresh, Sneha
AU - Romeo, Daniele
AU - Schunter, Celia
TI - Molecular Drivers of Mutualistic Association Between Anemone and Anemonefish
T2 - Genome biology and evolution
J2 - Genome Biol Evol
PY - 2026
DA - 2026/08/01
VL - 18
IS - 8
SP - evag197
SN - 1759-6653
PB - Oxford University Press
DO - 10.1093/gbe/evag197
UR - https://doi.org/10.1093/gbe/evag197
LA - en
ER -

CSL-JSON

{
"id": "10.1093/gbe/evag197",
"type": "article-journal",
"title": "Molecular Drivers of Mutualistic Association Between Anemone and Anemonefish",
"container-title": "Genome biology and evolution",
"author": [
{
"family": "Suresh",
"given": "Sneha"
},
{
"family": "Romeo",
"given": "Daniele"
},
{
"family": "Schunter",
"given": "Celia"
}
],
"container-title-short": "Genome Biol Evol",
"volume": "18",
"issue": "8",
"page": "evag197",
"DOI": "10.1093/gbe/evag197",
"PMID": "42535633",
"PMCID": "PMC13481981",
"ISSN": "1759-6653",
"publisher": "Oxford University Press",
"URL": "https://doi.org/10.1093/gbe/evag197",
"language": "en",
"issued": {
"date-parts": [
[
2026,
8,
1
]
]
}
}

The tracing map gets a citation of its own once an author has validated it and it has a DOI.

Similar papers

The papers with a page that share the most with this one: the tools found in their code, their categories, datasets, cited references and authors, the rarest counting most.

[1] doi:10.1093/nar/gkag294 [code]
Developmental stage dominates cell-type identity and reveals a chromatin regulatory function for Rad50 in Drosophila.
Journal: Nucleic acids research
In common: DESeq2, pheatmap, ggpubr, 2 other tools, 4 references
[2] doi:10.1038/s41467-026-76688-w [code]
Transcriptome profiling of human hypothalamic agouti-related protein and proopiomelanocortin neurons regulating energy homeostasis.
Journal: Nature communications
In common: DESeq2, pheatmap, ggpubr, 2 other tools, genetics / omics, cellular / molecular, 3 references
[3] doi:10.1038/s41586-026-10629-x [code]
Whole-genome duplication shaped cell-type evolution in the vertebrate brain.
Journal: Nature
In common: DESeq2, pheatmap, ggpubr, 2 other tools, other, genetics / omics, cellular / molecular, 2 references
[4] doi:10.7554/elife.107393 [code]
Chromosome-scale genome assembly of the European common cuttlefish &lt;i&gt;Sepia officinalis&lt;/i&gt;.
Journal: eLife
In common: DESeq2, pheatmap, ggplot2, 1 other tool, other, cellular / molecular, 3 references
[5] doi:10.1093/nar/gkag788 [code]
Neuronal activity-driven 3D chromatin dynamics in cortical pyramidal neurons depend on SATB2.
Journal: Nucleic acids research
In common: DESeq2, pheatmap, ggplot2, 1 other tool, cellular / molecular, 4 references
[6] doi:10.1038/s42003-026-10402-w [code]
Temporal orchestration of transcriptional and epigenomic programming underlying maternal embryonic diapause in a cricket model.
Journal: Communications biology
In common: DESeq2, pheatmap, ggplot2, 1 other tool, other, genetics / omics, 3 references
[7] doi:10.1016/j.stemcr.2026.102930 [code]
ZFHX4 is necessary for dopaminergic neuron differentiation and controls cell cycle by regulating LIN28A.
Journal: Stem cell reports
In common: DESeq2, pheatmap, ggpubr, 2 other tools, genetics / omics, cellular / molecular, 2 references
[8] doi:10.1038/s41586-026-10512-9 [code]
Astrocyte glucocorticoid receptor signalling restricts neuronal plasticity.
Journal: Nature
In common: DESeq2, pheatmap, ggpubr, 2 other tools, cellular / molecular, 2 references
[9] doi:10.21203/rs.3.rs-9927928/v1 [code]
Genome-wide and allele-resolved maps of the radial architecture of the mouse genome
Journal: Research Square (preprint)
In common: DESeq2, pheatmap, ggpubr, 2 other tools, genetics / omics, 2 references
[10] doi:10.1093/bioinformatics/btag592 [code]
Network-based stratification of allele-specific expression reveals patient subgroups in Huntington's disease.
Journal: Bioinformatics (Oxford, England)
In common: DESeq2, pheatmap, ggpubr, 2 other tools, genetics / omics, 2 references

Contribute

The authors of this paper can claim it, correct its record and validate its tracing map, and the maintainers of its code (its owner, or a public member of its organization) correct what it says of their repository; anyone signed in can ask for its removal. Every request goes to OSCR's own machine, which answers it; your account page follows them.

Sign in with ORCID to claim this paper as one of its authors, correct its record or validate its tracing map: when the paper's metadata lists your ORCID iD, you are recognized at once. Maintainers of its code: sign in with GitHub, then claim the repository on your account page.

Request its removal

To ask OSCR to remove this record, the copies of its authors' scripts or its tracing map, use the removal request page: signed in, you say who you are, what to remove and why, then review and confirm the request. Published rules decide every request (how).

Discussion, reproductions, activity

Discussion: questions and error reports about this paper and its code, from signed-in readers and its authors. It opens with sign-in.

Reproductions: reports from readers who ran the authors' code: what they reproduced, with which environment, commit and data. It opens with sign-in.

Activity: what happens around this paper: new versions of its record, its map's validation, discussions and reproductions. It opens with sign-in.