Preparation

Load packages and set options

Setup the default knitr options for the R code sections below.

Load the packages used throughout this analysis.

setwd("/share/users/Mike/AutismPaper/RNASeq/R")
suppressMessages({
  library('kableExtra')
  library('formattable')
  library('goseq')
  library('readr')
  library('stringr')
  library('tibble')
  library('dplyr')
  library('tidyr')
  library('tximport')
  library('DESeq2')
  library('ggplot2')
  library('BiocParallel')
  library('RColorBrewer')
  library('colorspace')
  # library('gplots')
  # library('superheat')
  library('corrplot')
  # library('pheatmap')
  library('VennDiagram')
  library('ggrepel')
  # library('gridExtra')
  # library('BatchQC')
  library('scales')
  register(MulticoreParam(detectCores()))
  register(SerialParam())
  options(tibble.print_max=100)
})

Import metadata

This section assembles the metadata for the datasets included in this analysis.

meta_int <- read_tsv('../../metadata/RNASeq.internal.tsv') %>%
  mutate(library = 'paired end') %>%
  mutate(readlength = 101) %>%
  mutate(tissue = 'SVZ') %>%
  mutate(path = paste0('/share/users/Mike/AutismPaper/RNASeq/OurData/', 
                           sample, '.refseqnorrna.salmon.paired/quant.sf')
  )

meta_ext <- read_tsv('../../metadata/RNASeq.external.tsv') %>%
  mutate(library = 'single end') %>%
  mutate(readlength = 75) %>%
  mutate(condition = 'Control') %>%
  mutate(path = paste0('/share/users/Mike/AutismPaper/RNASeq/external/', 
                           sample, '.refseqnorrna.salmon.single/quant.sf')
  )

meta_ext1 <-read_tsv('../../metadata/RNASeq.external2.tsv') %>%
  mutate(path = paste0('/share/users/Mike/AutismPaper/RNASeq/external2/', 
                           sample, '.refseqnorrna.salmon.single/quant.sf') 
)

metadata <- bind_rows(meta_int, meta_ext, meta_ext1) %>%
  mutate(age = factor(if_else(!is.na(age_wpc), 'Fetal', 
                              if_else(between(age_years,1,19),'Child','Adult')),
                      levels = c('Fetal', 'Child', 'Adult')
                      )
         ) %>%
  mutate(tissue1 = if_else(grepl('SVZ$', tissue), 'SVZ', 
                           if_else(grepl('frontal', tissue), 'FC',
                                   if_else(grepl('temporal', tissue), 'TC',
                                           tissue)
                                   )
                           )
         ) %>%
  mutate(sample = if_else(grepl('OurData', path),
                          paste(sep="_",sub("(\\w).*", "\\1", condition), 
                                paste0(age_years, 'y'),
                                sample), 
                          if_else(grepl('external2',path), 
                                  paste(sep='_', sub("(\\w).*", "\\1", condition), 
                                        tissue1, 
                                        paste0(age_years, 'y'), 
                                        sub(".*(\\d{2})$", "\\1", sample)
                                        ),
                                  paste(sep='_', tissue, 
                                        paste0(age_wpc, 'w'), 
                                        sub(".*(\\d{2})$", "\\1", sample)
                                        )
                                  )
                          )
  ) %>%
  filter(is.na(sex) | sex == 'M')

SVZ tissue analysis

DESeq2 analysis

Select SVZ metadata

Select the metadata corresponding to the SVZ tissue samples:

metadata1 <- metadata %>%
  filter(grepl('SVZ',tissue) & !grepl('H16$', sample)) %>%
  mutate(condition = factor(condition, levels=c('Control', 'Autism'))) %>%
  mutate(grp = factor(paste0(age,condition), levels=c('FetalControl', 'ChildControl', 'ChildAutism', 'AdultControl', 'AdultAutism'))) %>%
  arrange(sample) %>%
  dplyr::select( c('sample', 'condition', 'library', 'readlength', 'tissue', 'age', 'grp', 'path'))
metadata1

Read RefSeq transcript-to-genes data

Read the transcript-to-genes metadata table that associates RefSeq accessions (e.g. NM_xxxxxx) to Entrez Gene IDs, symbols, etc.

ttgfile <- '/share/db/refseq/human/refseq.108.rna.t2g.txt'
ttg <- read_tsv(ttgfile)

Read salmon output files

Read the files with the salmon estimated counts for each of the transcripts in RefSeq version 108 (http://www.ncbi.nlm.nih.gov/genome/annotation_euk/Homo_sapiens/108/).

files <- metadata1 %>% dplyr::select('sample', 'path') %>% spread('sample', 'path') %>% as.list() %>% unlist()
txi.salmon <- tximport(files, type = "salmon", tx2gene = ttg, txOut = TRUE, varReduce = TRUE)
metadata1.df <- metadata1 %>%  column_to_rownames(var="sample") %>% as.data.frame()
#rownames(txi.salmon$counts) <- as.character(rownames(txi.salmon$counts))

Perform Wald test

Create a DESeq dataset from the imported estimated counts object with a model design based on a factor made by combining the age group and the condition (grp) with a baseline level corresponding to the Control condition and Fetal age group.
Filter out all the rows that have at least 1 estimated read in the fetal samples, at least 1 in the control samples and at least 1 in the autism samples .
Perform a Wald test on this model.
List the available contrasts

dds <- DESeqDataSetFromTximport(txi = txi.salmon, colData = metadata1.df, design = ~grp)
c <- counts(dds) %>% as.tibble()
c[['target_id']] <- rownames(counts(dds))
c1 <- c %>% 
  gather(colnames(dds), key = 'sample', value = 'est_counts') %>% 
  left_join(metadata1 %>% select(c('sample', 'grp'))) %>% 
  dplyr::count(target_id, grp)
c2 <- c %>% 
  gather(colnames(dds), key = 'sample', value = 'est_counts') %>% 
  left_join(metadata1 %>% select(c('sample', 'grp'))) %>% 
  group_by(target_id, grp) %>% 
  summarise(lowvals =sum(if_else(est_counts < 5, 1, 0)))
c3 <- c1 %>% 
  inner_join(c2, by = c('target_id', 'grp')) %>% 
  mutate(notdetect = if_else(lowvals > n / 2, TRUE, FALSE))
c4 <- c3 %>% 
  dplyr::select(c('target_id', 'grp', 'notdetect'))  %>% 
  spread(key = 'grp', value = 'notdetect') %>%
  filter_at(vars(levels(metadata1$grp)), all_vars(!.))
dds <- dds[c4$target_id,]
rm(list=c('c','c1','c2','c3','c4'))
ddsWald <- DESeq(dds, test = "Wald", betaPrior = FALSE, parallel = TRUE)
resultsNames(ddsWald)

## [1] "Intercept"                        "grp_ChildControl_vs_FetalControl"
## [3] "grp_ChildAutism_vs_FetalControl"  "grp_AdultControl_vs_FetalControl"
## [5] "grp_AdultAutism_vs_FetalControl"

DESeq results

Extract the results from the model.

resultsChild   <- results(ddsWald, alpha = 0.05, lfcThreshold = 1, 
                          contrast = c('grp', 'ChildControl', 'FetalControl'))
resultsChild   <- lfcShrink(dds = ddsWald, 
                            res = resultsChild, 
                            contrast = c('grp', 'ChildControl', 'FetalControl'))

resultsAdult   <- results(ddsWald, alpha = 0.05, lfcThreshold = 1, 
                          contrast = c('grp', 'AdultControl', 'FetalControl'))
resultsAdult   <- lfcShrink(dds = ddsWald, 
                            res = resultsAdult, 
                            contrast = c('grp', 'AdultControl', 'FetalControl'))

resultsChildA  <- results(ddsWald, alpha = 0.05, lfcThreshold = 1, 
                          contrast = c('grp', 'ChildAutism', 'FetalControl'))
resultsChildA   <- lfcShrink(dds = ddsWald, 
                            res = resultsChildA, 
                            contrast = c('grp', 'ChildAutism', 'FetalControl'))

resultsChildArev  <- results(ddsWald, alpha = 0.05, lfcThreshold = 1, 
                          contrast = c('grp', 'ChildAutism', 'FetalControl'),
                          altHypothesis = 'lessAbs')
resultsChildArev   <- lfcShrink(dds = ddsWald, 
                                res = resultsChildArev, 
                                contrast = c('grp', 'ChildAutism', 'FetalControl'))

resultsAdultA  <- results(ddsWald, alpha = 0.05, lfcThreshold = 1, 
                          contrast = c('grp', 'AdultAutism', 'FetalControl'))
resultsAdultA   <- lfcShrink(dds = ddsWald, 
                            res = resultsAdultA, 
                            contrast = c('grp', 'AdultAutism', 'FetalControl'))

signifChild <- resultsChild %>%
  as.tibble() %>% rownames_to_column(var='target_id') %>% 
  filter(padj<0.05) %>% arrange(padj) %>% 
  left_join(ttg, by = 'target_id') %>% 
  dplyr::select(c('target_id', 'baseMean', 'log2FoldChange', 'lfcSE', 'pvalue', 'padj', 'gene_id', 'symbol', 'description', 'status', 'molecule_type'))

signifAdult <- resultsAdult %>%
  as.tibble(.) %>% rownames_to_column(var='target_id') %>%
  filter(padj<0.05) %>% arrange(padj) %>%
  left_join(ttg, by = 'target_id') %>%
  dplyr::select(c('target_id', 'baseMean', 'log2FoldChange', 'lfcSE', 'pvalue', 'padj', 'gene_id', 'symbol', 'description', 'status', 'molecule_type'))

signifChildA <- resultsChildA %>%
  as.tibble(.) %>% rownames_to_column(var='target_id') %>%
  filter(padj<0.05) %>% arrange(padj) %>%
  left_join(ttg, by = 'target_id') %>%
  dplyr::select(c('target_id', 'baseMean', 'log2FoldChange', 'lfcSE', 'pvalue', 'padj', 'gene_id', 'symbol', 'description', 'status', 'molecule_type'))

signifChildArev <- resultsChildArev %>%
  as.tibble(.) %>% rownames_to_column(var='target_id') %>%
  filter(padj<0.05) %>% arrange(padj) %>%
  left_join(ttg, by = 'target_id') %>%
  dplyr::select(c('target_id', 'baseMean', 'log2FoldChange', 'lfcSE', 'pvalue', 'padj', 'gene_id', 'symbol', 'description', 'status', 'molecule_type'))

signifAdultA <- resultsAdultA %>%
  as.tibble(.) %>% rownames_to_column(var='target_id') %>%
  filter(padj<0.05) %>% arrange(padj) %>%
  left_join(ttg, by = 'target_id') %>%
  dplyr::select(c('target_id', 'baseMean', 'log2FoldChange', 'lfcSE', 'pvalue', 'padj', 'gene_id', 'symbol', 'description', 'status', 'molecule_type'))

Control - Fetal

Summary of the number of transcripts and genes differentially expressed in control samples of different age groups vs. the fetal samples (“developmental genes”):

rsum_control_fetal <- tribble( ~Counts, ~Child, ~Adult,
                               'Transcripts', 
                               length(signifChild$target_id), 
                               length(signifAdult$target_id), 
                               'Genes', 
                               length(unique(signifChild$gene_id)),
                               length(unique(signifAdult$gene_id))
                               )
kable(rsum_control_fetal, format = 'markdown')

Counts	Child	Adult
Transcripts	1942	2361
Genes	1770	2092

Venn diagram:

grid.draw(venn.diagram(list(Child=signifChild$target_id, Adult=signifAdult$target_id), 
             filename = NULL,
             alpha = 0.5,
             inverted = TRUE, 
             cat.pos = c(45, 315),
             fill = c("cornflowerblue", "darkorchid1"),
             cat.col = c("darkblue", "darkorchid4")
             )
)

Autism - Fetal

Summary of the number of transcripts and genes differentially expressed in autism samples of different age groups vs. the fetal samples:

rsum_autism_fetal <- tribble( ~Counts, ~Child, ~Adult, 
                               'Transcripts', 
                               length(signifChildA$target_id), 
                               length(signifAdultA$target_id), 
                               'Genes', 
                               length(unique(signifChildA$gene_id)),
                               length(unique(signifAdultA$gene_id))
                               )
kable(rsum_autism_fetal, format = 'markdown')

Counts	Child	Adult
Transcripts	2653	2118
Genes	2384	1904

Venn diagram:

grid.draw(venn.diagram(list(Child=signifChildA$target_id, Adult=signifAdultA$target_id),
             filename = NULL,
             alpha = 0.5,
             fill = c("cornflowerblue", "darkorchid1"),
             cat.col = c("darkblue", "darkorchid4")
             )
)

QC

PCA

## Transform count data
intgroup = "grp"
rld <- rlog(ddsWald, blind = FALSE)

## Perform PCA analysis and get percent of variance explained
data_pca <- plotPCA(rld, intgroup = intgroup, ntop = 1000, returnData = TRUE)
percentVar <- round(100 * attr(data_pca, "percentVar"))

## Make plot
data_pca %>%
ggplot(aes_string(x = "PC1", y = "PC2", color = "grp")) +
 geom_point(size = 3) +
 xlab(paste0("PC1: ", percentVar[1], "% variance")) +
 ylab(paste0("PC2: ", percentVar[2], "% variance")) +
 coord_fixed() +
 geom_text(aes(label=name), nudge_x = 1, hjust = 0, size = 3, check_overlap = FALSE)

The above plot shows the first two principal components that explain the variability in the data using the regularized log count data. If you are unfamiliar with principal component analysis, you might want to check the Wikipedia entry or this interactive explanation. In this case, the first and second principal component explain percentVar[1] and percentVar[2] percent of the variance respectively.

Volcano plots

Control - Fetal

ChildControl - Fetal

resultsChild %>%
  as.tibble()  %>%
  rownames_to_column(var='target_id') %>%
  mutate(sig=ifelse(target_id %in% signifChild[['target_id']], "FDR<0.05", "Not Sig")) %>%
  left_join(ttg, by = 'target_id') %>% 
  ggplot(aes(log2FoldChange, -log10(padj))) + 
  geom_point(aes(col=sig), alpha=1/10) + 
  scale_color_manual(values=c('red', 'black'))

AdultControl - Fetal

resultsAdult %>%
  as.tibble()  %>%
  rownames_to_column(var='target_id') %>%
  mutate(sig=ifelse(target_id %in% signifAdult[['target_id']], "FDR<0.05", "Not Sig")) %>%
  left_join(ttg, by = 'target_id') %>% 
  ggplot(aes(log2FoldChange, -log10(padj))) + 
  geom_point(aes(col=sig), alpha=1/10) + 
  scale_color_manual(values=c("red", "black"))

Autism - Fetal

ChildAutism - Fetal

resultsChildA %>%
  as.tibble()  %>%
  rownames_to_column(var='target_id') %>%
  mutate(sig=ifelse(target_id %in% signifChildA[['target_id']], "FDR<0.05", "Not Sig")) %>%
  left_join(ttg, by = 'target_id') %>% 
  ggplot(aes(log2FoldChange, -log10(padj))) + 
  geom_point(aes(col=sig), alpha=1/10) + 
  scale_color_manual(values=c("red", "black"))

AdultAutism - Fetal

resultsAdultA %>%
  as.tibble()  %>%
  rownames_to_column(var='target_id') %>%
  mutate(sig=ifelse(target_id %in% signifAdultA[['target_id']], "FDR<0.05", "Not Sig")) %>%
  left_join(ttg, by = 'target_id') %>% 
  ggplot(aes(log2FoldChange, -log10(padj))) + 
  geom_point(aes(col=sig), alpha=1/10) + 
  scale_color_manual(values=c("red", "black"))

Gene lists

Control - Fetal

ChildControl - Fetal

This table shows the 50 most significant genes (out of 1770) differentially expressed between control child samples and fetal samples.

genes_control_fetal_child <- signifChild %>%
  mutate(gene = text_spec(symbol, link = paste0('http://www.genecards.org/cgi-bin/carddisp.pl?gene=', symbol)),
         log2FoldChange = ifelse(log2FoldChange < 0,
                                 color_tile("red", "white")(signif(log2FoldChange, 3)),
                                 color_tile("white", "darkgreen")(signif(log2FoldChange, 3))),
         padj = format(padj, digits=3, scientific = TRUE),
         transcript = target_id)
genes_control_fetal_child  %>%
  dplyr::select(c('gene', 'transcript', 'log2FoldChange', 'padj', 'description', 'status', 'molecule_type')) %>%
  head(n = 50) %>%
  kable("html", escape = F) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),  font_size = 8, full_width = F)

gene	transcript	log2FoldChange	padj	description	status	molecule_type
SOX11	NM_003108.3	-7.34	1.63e-218	SRY-box 11 (SOX11)	REVIEWED	mRNA
MEX3A	NM_001093725.1	-6.01	1.55e-108	mex-3 RNA binding family member A (MEX3A)	VALIDATED	mRNA
GPC2	NM_152742.2	-5.30	2.29e-107	glypican 2 (GPC2)	VALIDATED	mRNA
KIF11	NM_004523.3	-5.25	2.42e-83	kinesin family member 11 (KIF11)	REVIEWED	mRNA
SOX4	NM_003107.2	-5.54	8.93e-81	SRY-box 4 (SOX4)	REVIEWED	mRNA
IGIP	NM_001007189.1	3.63	2.16e-62	IgA inducing protein (IGIP)	VALIDATED	mRNA
CDK1	NM_001786.4	-7.95	5.15e-61	cyclin dependent kinase 1 (CDK1), transcript variant 1	REVIEWED	mRNA
KNL1	XM_017022432.1	-6.00	3.57e-58	PREDICTED: Homo sapiens cancer susceptibility candidate 5 (CASC5), transcript variant X1	MODEL	mRNA
SPC25	NM_020675.3	-5.18	5.93e-57	SPC25, NDC80 kinetochore complex component (SPC25)	REVIEWED	mRNA
TMEM151A	NM_153266.3	7.28	3.68e-52	transmembrane protein 151A (TMEM151A)	VALIDATED	mRNA
TPPP	NM_007030.2	6.73	1.30e-51	tubulin polymerization promoting protein (TPPP)	VALIDATED	mRNA
RRM2	NM_001034.3	-6.55	5.56e-50	ribonucleotide reductase regulatory subunit M2 (RRM2), transcript variant 2	REVIEWED	mRNA
SNORA109	NR_132964.1	6.70	3.19e-48	small nucleolar RNA, H/ACA box 109 (SNORA109)	VALIDATED	small nucleolar RNA
CCNB2	NM_004701.3	-7.50	8.26e-46	cyclin B2 (CCNB2)	REVIEWED	mRNA
PLEKHB1	NM_001130035.1	6.35	4.31e-45	pleckstrin homology domain containing B1 (PLEKHB1), transcript variant 4	VALIDATED	mRNA
HCN2	NM_001194.3	7.82	7.41e-45	hyperpolarization activated cyclic nucleotide gated potassium and sodium channel 2 (HCN2)	REVIEWED	mRNA
ZBTB47	NM_145166.3	3.81	1.13e-43	zinc finger and BTB domain containing 47 (ZBTB47)	VALIDATED	mRNA
GDF1	NM_001492.5	4.60	4.32e-43	growth differentiation factor 1 (GDF1)	REVIEWED	mRNA
CERS1	NM_021267.4	4.60	4.32e-43	ceramide synthase 1 (CERS1), transcript variant 1	REVIEWED	mRNA
LMNB1	NM_005573.3	-5.07	4.00e-41	lamin B1 (LMNB1), transcript variant 1	REVIEWED	mRNA
PTTG1	NM_004219.3	-5.58	1.68e-40	pituitary tumor-transforming 1 (PTTG1), transcript variant 2	REVIEWED	mRNA
ZNF365	NM_014951.2	6.14	7.76e-38	zinc finger protein 365 (ZNF365), transcript variant A	REVIEWED	mRNA
CCDC85B	NM_006848.2	5.23	6.85e-34	coiled-coil domain containing 85B (CCDC85B)	VALIDATED	mRNA
TF	NM_001063.3	9.21	2.93e-33	transferrin (TF)	REVIEWED	mRNA
ENDOD1	NM_015036.2	5.59	1.03e-32	endonuclease domain containing 1 (ENDOD1)	VALIDATED	mRNA
HMGB3	NM_005342.3	-3.59	1.10e-31	high mobility group box 3 (HMGB3), transcript variant 2	REVIEWED	mRNA
PLLP	NM_015993.2	7.25	2.58e-31	plasmolipin (PLLP)	VALIDATED	mRNA
EPAS1	NM_001430.4	3.76	6.03e-31	endothelial PAS domain protein 1 (EPAS1)	REVIEWED	mRNA
HSPA2	NM_021979.3	6.61	1.78e-30	heat shock protein family A (Hsp70) member 2 (HSPA2)	VALIDATED	mRNA
TIAM2	NM_012454.3	-6.26	3.90e-30	T-cell lymphoma invasion and metastasis 2 (TIAM2), transcript variant 1	REVIEWED	mRNA
LINC01158	NR_131235.1	-5.88	6.27e-30	long intergenic non-protein coding RNA 1158 (LINC01158), transcript variant 4	VALIDATED	long non-coding RNA
TKTL1	NM_001145934.1	-6.03	2.31e-29	transketolase like 1 (TKTL1), transcript variant 3	REVIEWED	mRNA
TMPO	NM_001032284.2	-3.75	5.28e-29	thymopoietin (TMPO), transcript variant 3	REVIEWED	mRNA
PINK1	NM_032409.2	3.17	1.81e-28	PTEN induced putative kinase 1 (PINK1)	REVIEWED	mRNA
CST3	NM_000099.3	4.93	4.08e-28	cystatin C (CST3), transcript variant 1	REVIEWED	mRNA
PLXNB3	NM_005393.2	6.47	8.42e-28	plexin B3 (PLXNB3), transcript variant 1	REVIEWED	mRNA
GINS2	NM_016095.2	-4.59	1.68e-27	GINS complex subunit 2 (GINS2)	VALIDATED	mRNA
FABP7	NM_001319041.1	-5.11	1.73e-27	fatty acid binding protein 7 (FABP7), transcript variant 3	REVIEWED	mRNA
FAIM2	NM_012306.3	5.50	3.69e-27	Fas apoptotic inhibitory molecule 2 (FAIM2)	VALIDATED	mRNA
NAP1L1	XR_001748714.1	-3.77	4.78e-27	PREDICTED: Homo sapiens nucleosome assembly protein 1 like 1 (NAP1L1), transcript variant X4	MODEL	misc_RNA
MAD2L1	NM_002358.3	-4.01	1.41e-26	mitotic arrest deficient 2 like 1 (MAD2L1)	REVIEWED	mRNA
CNP	NM_033133.4	5.66	1.77e-26	2’,3’-cyclic nucleotide 3’ phosphodiesterase (CNP), transcript variant 1	VALIDATED	mRNA
NEUROD2	NM_006160.3	-4.61	2.50e-26	neuronal differentiation 2 (NEUROD2)	REVIEWED	mRNA
PRNP	NM_001080123.2	3.87	5.87e-26	prion protein (PRNP), transcript variant 5	REVIEWED	mRNA
SPNS2	NM_001124758.2	4.70	5.87e-26	sphingolipid transporter 2 (SPNS2)	REVIEWED	mRNA
QSER1	NM_001076786.2	-3.07	1.03e-25	glutamine and serine rich 1 (QSER1)	VALIDATED	mRNA
FOXO3B	NR_026718.1	-3.11	1.57e-25	forkhead box O3B pseudogene (FOXO3B)	PROVISIONAL	non-coding RNA
TMEM161B-AS1	NR_105020.1	-5.16	1.57e-25	TMEM161B antisense RNA 1 (TMEM161B-AS1), transcript variant 5	VALIDATED	long non-coding RNA
MT3	NM_005954.2	4.97	2.21e-25	metallothionein 3 (MT3)	VALIDATED	mRNA
MCM2	NM_004526.3	-4.19	2.37e-25	minichromosome maintenance complex component 2 (MCM2), transcript variant 1	REVIEWED	mRNA

AdultControl - Fetal

This table shows the 50 most significant genes (out of 2092) differentially expressed between control adult samples and fetal samples.

genes_control_fetal_adult <- signifAdult %>%
  mutate(gene = text_spec(symbol, link = paste0('http://www.genecards.org/cgi-bin/carddisp.pl?gene=', symbol)),
         log2FoldChange = ifelse(log2FoldChange < 0,
                                 color_tile("red", "white")(signif(log2FoldChange, 3)),
                                 color_tile("white", "darkgreen")(signif(log2FoldChange, 3))),
         padj = format(padj, digits=3, scientific = TRUE),
         transcript = target_id)
genes_control_fetal_adult %>%
  dplyr::select(c('gene', 'transcript', 'log2FoldChange', 'padj', 'description', 'status', 'molecule_type')) %>%
  head(n = 50) %>%
  kable("html", escape = F) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),  font_size = 8, full_width = F)

gene	transcript	log2FoldChange	padj	description	status	molecule_type
SOX11	NM_003108.3	-8.170	1.83e-262	SRY-box 11 (SOX11)	REVIEWED	mRNA
SOX4	NM_003107.2	-6.550	5.00e-117	SRY-box 4 (SOX4)	REVIEWED	mRNA
MEX3A	NM_001093725.1	-6.340	1.26e-114	mex-3 RNA binding family member A (MEX3A)	VALIDATED	mRNA
KIF11	NM_004523.3	-5.260	3.20e-83	kinesin family member 11 (KIF11)	REVIEWED	mRNA
GPC2	NM_152742.2	-4.430	1.69e-76	glypican 2 (GPC2)	VALIDATED	mRNA
IGIP	NM_001007189.1	3.870	4.26e-75	IgA inducing protein (IGIP)	VALIDATED	mRNA
TMEM151A	NM_153266.3	8.000	3.76e-65	transmembrane protein 151A (TMEM151A)	VALIDATED	mRNA
TPPP	NM_007030.2	7.380	1.04e-64	tubulin polymerization promoting protein (TPPP)	VALIDATED	mRNA
CDK1	NM_001786.4	-7.490	6.44e-61	cyclin dependent kinase 1 (CDK1), transcript variant 1	REVIEWED	mRNA
PLEKHB1	NM_001130035.1	6.970	1.08e-56	pleckstrin homology domain containing B1 (PLEKHB1), transcript variant 4	VALIDATED	mRNA
KNL1	XM_017022432.1	-5.810	5.58e-56	PREDICTED: Homo sapiens cancer susceptibility candidate 5 (CASC5), transcript variant X1	MODEL	mRNA
RRM2	NM_001034.3	-7.270	1.63e-55	ribonucleotide reductase regulatory subunit M2 (RRM2), transcript variant 2	REVIEWED	mRNA
ENDOD1	NM_015036.2	6.770	1.19e-52	endonuclease domain containing 1 (ENDOD1)	VALIDATED	mRNA
HMGB3	NM_005342.3	-4.340	4.26e-51	high mobility group box 3 (HMGB3), transcript variant 2	REVIEWED	mRNA
LMNB1	NM_005573.3	-5.750	5.70e-50	lamin B1 (LMNB1), transcript variant 1	REVIEWED	mRNA
TF	NM_001063.3	10.600	8.16e-47	transferrin (TF)	REVIEWED	mRNA
SPC25	NM_020675.3	-4.370	1.21e-45	SPC25, NDC80 kinetochore complex component (SPC25)	REVIEWED	mRNA
PTGFRN	XM_017001874.1	-3.680	1.99e-44	PREDICTED: Homo sapiens prostaglandin F2 receptor inhibitor (PTGFRN), transcript variant X1	MODEL	mRNA
HSPA2	NM_021979.3	7.690	3.70e-44	heat shock protein family A (Hsp70) member 2 (HSPA2)	VALIDATED	mRNA
PLLP	NM_015993.2	8.350	5.97e-44	plasmolipin (PLLP)	VALIDATED	mRNA
ZBTB47	NM_145166.3	3.790	3.34e-43	zinc finger and BTB domain containing 47 (ZBTB47)	VALIDATED	mRNA
GDF1	NM_001492.5	4.610	3.78e-43	growth differentiation factor 1 (GDF1)	REVIEWED	mRNA
CERS1	NM_021267.4	4.610	3.78e-43	ceramide synthase 1 (CERS1), transcript variant 1	REVIEWED	mRNA
PTTG1	NM_004219.3	-6.160	1.35e-42	pituitary tumor-transforming 1 (PTTG1), transcript variant 2	REVIEWED	mRNA
PLXNB3	NM_005393.2	7.690	1.80e-42	plexin B3 (PLXNB3), transcript variant 1	REVIEWED	mRNA
HCN2	NM_001194.3	7.580	7.22e-42	hyperpolarization activated cyclic nucleotide gated potassium and sodium channel 2 (HCN2)	REVIEWED	mRNA
ZNF365	NM_014951.2	6.370	1.24e-41	zinc finger protein 365 (ZNF365), transcript variant A	REVIEWED	mRNA
DCHS1	NM_003737.3	-3.330	2.66e-41	dachsous cadherin-related 1 (DCHS1)	REVIEWED	mRNA
SPNS2	NM_001124758.2	5.490	7.24e-39	sphingolipid transporter 2 (SPNS2)	REVIEWED	mRNA
CNP	NM_033133.4	6.550	2.85e-38	2’,3’-cyclic nucleotide 3’ phosphodiesterase (CNP), transcript variant 1	VALIDATED	mRNA
SNORA109	NR_132964.1	6.070	5.53e-38	small nucleolar RNA, H/ACA box 109 (SNORA109)	VALIDATED	small nucleolar RNA
APLP1	NM_005166.4	4.630	9.78e-37	amyloid beta precursor like protein 1 (APLP1), transcript variant 2	VALIDATED	mRNA
TIAM2	NM_012454.3	-7.050	2.47e-35	T-cell lymphoma invasion and metastasis 2 (TIAM2), transcript variant 1	REVIEWED	mRNA
ENPP4	NM_014936.4	5.120	5.39e-35	ectonucleotide pyrophosphatase/phosphodiesterase 4 (putative) (ENPP4)	VALIDATED	mRNA
CCNB2	NM_004701.3	-5.350	6.93e-35	cyclin B2 (CCNB2)	REVIEWED	mRNA
GPR37	NM_005302.3	7.830	1.08e-34	G protein-coupled receptor 37 (GPR37)	REVIEWED	mRNA
LINC01158	NR_131235.1	-6.270	6.26e-34	long intergenic non-protein coding RNA 1158 (LINC01158), transcript variant 4	VALIDATED	long non-coding RNA
SPTBN1	NM_178313.2	5.870	2.71e-32	spectrin beta, non-erythrocytic 1 (SPTBN1), transcript variant 2	REVIEWED	mRNA
CKS2	NM_001827.2	-5.730	3.12e-32	CDC28 protein kinase regulatory subunit 2 (CKS2)	REVIEWED	mRNA
TMEM161B-AS1	NR_105020.1	-6.010	5.48e-32	TMEM161B antisense RNA 1 (TMEM161B-AS1), transcript variant 5	VALIDATED	long non-coding RNA
DBN1	NM_004395.3	-3.970	8.78e-32	drebrin 1 (DBN1), transcript variant 1	REVIEWED	mRNA
TACC3	NM_006342.2	-4.890	2.68e-31	transforming acidic coiled-coil containing protein 3 (TACC3)	REVIEWED	mRNA
EFNB2	NM_004093.3	-4.270	2.86e-31	ephrin B2 (EFNB2)	REVIEWED	mRNA
KIF1C	XM_005256424.2	5.620	5.34e-31	PREDICTED: Homo sapiens kinesin family member 1C (KIF1C), transcript variant X1	MODEL	mRNA
UBL3	NM_007106.3	2.520	5.92e-31	ubiquitin like 3 (UBL3)	VALIDATED	mRNA
MEX3B	NM_032246.4	-4.580	1.09e-29	mex-3 RNA binding family member B (MEX3B)	REVIEWED	mRNA
TMPO	NM_001032284.2	-3.730	1.31e-28	thymopoietin (TMPO), transcript variant 3	REVIEWED	mRNA
MAD2L1	NM_002358.3	-4.150	1.56e-28	mitotic arrest deficient 2 like 1 (MAD2L1)	REVIEWED	mRNA
FN3K	NM_022158.3	3.560	3.91e-28	fructosamine 3 kinase (FN3K)	VALIDATED	mRNA
TKTL1	NM_001145934.1	-6.380	1.40e-27	transketolase like 1 (TKTL1), transcript variant 3	REVIEWED	mRNA

Autism - Fetal

ChildAutism - Fetal

This table shows the 50 most significant genes (out of 2384) differentially expressed between autism child samples and fetal samples.

genes_autism_fetal_child <- signifChildA %>%
  mutate(gene = text_spec(symbol, link = paste0('http://www.genecards.org/cgi-bin/carddisp.pl?gene=', symbol)),
         log2FoldChange = ifelse(log2FoldChange < 0,
                                 color_tile("red", "white")(signif(log2FoldChange, 3)),
                                 color_tile("white", "darkgreen")(signif(log2FoldChange, 3))),
         padj = format(padj, digits=3, scientific = TRUE),
         transcript = target_id)
genes_autism_fetal_child %>%
  dplyr::select(c('gene', 'transcript', 'log2FoldChange', 'padj', 'description', 'status', 'molecule_type')) %>%
  head(n = 50) %>%
  kable("html", escape = F) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),  font_size = 8, full_width = F)

gene	transcript	log2FoldChange	padj	description	status	molecule_type
SOX11	NM_003108.3	-6.76	6.94e-229	SRY-box 11 (SOX11)	REVIEWED	mRNA
MEX3A	NM_001093725.1	-5.97	1.19e-129	mex-3 RNA binding family member A (MEX3A)	VALIDATED	mRNA
GPC2	NM_152742.2	-5.42	6.27e-127	glypican 2 (GPC2)	VALIDATED	mRNA
KIF11	NM_004523.3	-5.73	1.31e-113	kinesin family member 11 (KIF11)	REVIEWED	mRNA
SOX4	NM_003107.2	-5.39	2.04e-94	SRY-box 4 (SOX4)	REVIEWED	mRNA
IGIP	NM_001007189.1	3.65	1.35e-72	IgA inducing protein (IGIP)	VALIDATED	mRNA
SNORA109	NR_132964.1	7.15	6.43e-70	small nucleolar RNA, H/ACA box 109 (SNORA109)	VALIDATED	small nucleolar RNA
CDK1	NM_001786.4	-7.87	1.71e-69	cyclin dependent kinase 1 (CDK1), transcript variant 1	REVIEWED	mRNA
KNL1	XM_017022432.1	-6.42	1.55e-67	PREDICTED: Homo sapiens cancer susceptibility candidate 5 (CASC5), transcript variant X1	MODEL	mRNA
SPC25	NM_020675.3	-5.54	1.73e-65	SPC25, NDC80 kinetochore complex component (SPC25)	REVIEWED	mRNA
RRM2	NM_001034.3	-6.92	8.94e-64	ribonucleotide reductase regulatory subunit M2 (RRM2), transcript variant 2	REVIEWED	mRNA
HCN2	NM_001194.3	8.34	1.09e-62	hyperpolarization activated cyclic nucleotide gated potassium and sodium channel 2 (HCN2)	REVIEWED	mRNA
TMEM151A	NM_153266.3	7.27	2.63e-62	transmembrane protein 151A (TMEM151A)	VALIDATED	mRNA
GDF1	NM_001492.5	4.87	1.26e-60	growth differentiation factor 1 (GDF1)	REVIEWED	mRNA
CERS1	NM_021267.4	4.87	1.26e-60	ceramide synthase 1 (CERS1), transcript variant 1	REVIEWED	mRNA
TPPP	NM_007030.2	6.52	3.25e-59	tubulin polymerization promoting protein (TPPP)	VALIDATED	mRNA
LMNB1	NM_005573.3	-5.46	7.24e-57	lamin B1 (LMNB1), transcript variant 1	REVIEWED	mRNA
CCNB2	NM_004701.3	-6.79	5.80e-54	cyclin B2 (CCNB2)	REVIEWED	mRNA
PLEKHB1	NM_001130035.1	6.25	5.85e-54	pleckstrin homology domain containing B1 (PLEKHB1), transcript variant 4	VALIDATED	mRNA
ZBTB47	NM_145166.3	3.80	1.37e-53	zinc finger and BTB domain containing 47 (ZBTB47)	VALIDATED	mRNA
QSER1	NM_001076786.2	-3.60	3.95e-50	glutamine and serine rich 1 (QSER1)	VALIDATED	mRNA
ZNF365	NM_014951.2	6.38	3.88e-49	zinc finger protein 365 (ZNF365), transcript variant A	REVIEWED	mRNA
CCDC85B	NM_006848.2	5.52	1.77e-47	coiled-coil domain containing 85B (CCDC85B)	VALIDATED	mRNA
FOXO3B	NR_026718.1	-3.68	1.47e-46	forkhead box O3B pseudogene (FOXO3B)	PROVISIONAL	non-coding RNA
LINC01158	NR_131235.1	-6.53	3.31e-45	long intergenic non-protein coding RNA 1158 (LINC01158), transcript variant 4	VALIDATED	long non-coding RNA
DCHS1	NM_003737.3	-3.17	7.55e-45	dachsous cadherin-related 1 (DCHS1)	REVIEWED	mRNA
PINK1	NM_032409.2	3.40	7.48e-44	PTEN induced putative kinase 1 (PINK1)	REVIEWED	mRNA
TMPO	NM_001032284.2	-4.00	1.45e-42	thymopoietin (TMPO), transcript variant 3	REVIEWED	mRNA
NEUROD2	NM_006160.3	-5.16	4.31e-42	neuronal differentiation 2 (NEUROD2)	REVIEWED	mRNA
TF	NM_001063.3	9.28	1.77e-41	transferrin (TF)	REVIEWED	mRNA
FAIM2	NM_012306.3	5.82	5.50e-39	Fas apoptotic inhibitory molecule 2 (FAIM2)	VALIDATED	mRNA
MT3	NM_005954.2	5.39	8.32e-39	metallothionein 3 (MT3)	VALIDATED	mRNA
PTTG1	NM_004219.3	-4.85	1.23e-38	pituitary tumor-transforming 1 (PTTG1), transcript variant 2	REVIEWED	mRNA
TBR1	NM_006593.3	-8.09	3.19e-37	T-box, brain 1 (TBR1)	REVIEWED	mRNA
ATP1B1	NM_001677.3	4.26	1.52e-36	ATPase Na+/K+ transporting subunit beta 1 (ATP1B1)	REVIEWED	mRNA
ENDOD1	NM_015036.2	5.38	1.96e-36	endonuclease domain containing 1 (ENDOD1)	VALIDATED	mRNA
CST3	NM_000099.3	4.99	4.35e-36	cystatin C (CST3), transcript variant 1	REVIEWED	mRNA
HSPA2	NM_021979.3	6.49	5.03e-36	heat shock protein family A (Hsp70) member 2 (HSPA2)	VALIDATED	mRNA
TIAM2	NM_012454.3	-6.17	1.59e-35	T-cell lymphoma invasion and metastasis 2 (TIAM2), transcript variant 1	REVIEWED	mRNA
PLLP	NM_015993.2	7.06	3.19e-35	plasmolipin (PLLP)	VALIDATED	mRNA
HNRNPA1	NM_002136.3	-2.70	3.64e-35	heterogeneous nuclear ribonucleoprotein A1 (HNRNPA1), transcript variant 1	REVIEWED	mRNA
NAP1L1	XR_001748714.1	-3.80	1.30e-34	PREDICTED: Homo sapiens nucleosome assembly protein 1 like 1 (NAP1L1), transcript variant X4	MODEL	misc_RNA
SPRYD3	NM_032840.2	4.54	1.72e-33	SPRY domain containing 3 (SPRYD3)	VALIDATED	mRNA
TMEM161B-AS1	NR_105020.1	-5.38	2.11e-33	TMEM161B antisense RNA 1 (TMEM161B-AS1), transcript variant 5	VALIDATED	long non-coding RNA
GINS2	NM_016095.2	-4.62	2.58e-33	GINS complex subunit 2 (GINS2)	VALIDATED	mRNA
SETBP1	NM_015559.2	-3.82	3.55e-33	SET binding protein 1 (SETBP1), transcript variant 1	REVIEWED	mRNA
TKTL1	NM_001145934.1	-5.87	3.96e-33	transketolase like 1 (TKTL1), transcript variant 3	REVIEWED	mRNA
FABP7	NM_001319041.1	-5.07	9.21e-33	fatty acid binding protein 7 (FABP7), transcript variant 3	REVIEWED	mRNA
TACC3	NM_006342.2	-4.36	1.22e-32	transforming acidic coiled-coil containing protein 3 (TACC3)	REVIEWED	mRNA
HELLS	XR_001747095.1	-4.29	2.04e-32	PREDICTED: Homo sapiens helicase, lymphoid-specific (HELLS), transcript variant X1	MODEL	misc_RNA

AdultAutism - Fetal

This table shows the 50 most significant genes (out of 1904) differentially expressed between autism adult samples and fetal samples.

genes_autism_fetal_adult <- signifAdultA %>%
  mutate(gene = text_spec(symbol, link = paste0('http://www.genecards.org/cgi-bin/carddisp.pl?gene=', symbol)),
         log2FoldChange = ifelse(log2FoldChange < 0,
                                 color_tile("red", "white")(signif(log2FoldChange, 3)),
                                 color_tile("white", "darkgreen")(signif(log2FoldChange, 3))),
         padj = format(padj, digits=3, scientific = TRUE),
         transcript = target_id)
genes_autism_fetal_adult %>%
  dplyr::select(c('gene', 'transcript', 'log2FoldChange', 'padj', 'description', 'status', 'molecule_type')) %>%
  head(n = 50) %>%
  kable("html", escape = F) %>% 
  kable_styling(bootstrap_options = c("striped", "hover", "condensed"),  font_size = 8, full_width = F)

gene	transcript	log2FoldChange	padj	description	status	molecule_type
SOX11	NM_003108.3	-7.330	4.36e-218	SRY-box 11 (SOX11)	REVIEWED	mRNA
SOX4	NM_003107.2	-6.570	3.27e-118	SRY-box 4 (SOX4)	REVIEWED	mRNA
MEX3A	NM_001093725.1	-6.290	2.44e-113	mex-3 RNA binding family member A (MEX3A)	VALIDATED	mRNA
KIF11	NM_004523.3	-5.440	3.97e-88	kinesin family member 11 (KIF11)	REVIEWED	mRNA
GPC2	NM_152742.2	-4.810	3.97e-88	glypican 2 (GPC2)	VALIDATED	mRNA
IGIP	NM_001007189.1	3.840	3.73e-73	IgA inducing protein (IGIP)	VALIDATED	mRNA
TPPP	NM_007030.2	7.310	5.03e-63	tubulin polymerization promoting protein (TPPP)	VALIDATED	mRNA
TMEM151A	NM_153266.3	7.740	2.87e-60	transmembrane protein 151A (TMEM151A)	VALIDATED	mRNA
CDK1	NM_001786.4	-7.890	8.00e-60	cyclin dependent kinase 1 (CDK1), transcript variant 1	REVIEWED	mRNA
KNL1	XM_017022432.1	-5.700	2.12e-54	PREDICTED: Homo sapiens cancer susceptibility candidate 5 (CASC5), transcript variant X1	MODEL	mRNA
PLEKHB1	NM_001130035.1	6.760	1.04e-52	pleckstrin homology domain containing B1 (PLEKHB1), transcript variant 4	VALIDATED	mRNA
HMGB3	NM_005342.3	-4.350	3.32e-51	high mobility group box 3 (HMGB3), transcript variant 2	REVIEWED	mRNA
SPC25	NM_020675.3	-4.620	1.86e-49	SPC25, NDC80 kinetochore complex component (SPC25)	REVIEWED	mRNA
LMNB1	NM_005573.3	-5.620	2.14e-48	lamin B1 (LMNB1), transcript variant 1	REVIEWED	mRNA
RRM2	NM_001034.3	-6.400	5.97e-48	ribonucleotide reductase regulatory subunit M2 (RRM2), transcript variant 2	REVIEWED	mRNA
ZBTB47	NM_145166.3	3.800	1.99e-43	zinc finger and BTB domain containing 47 (ZBTB47)	VALIDATED	mRNA
PTTG1	NM_004219.3	-6.350	5.73e-43	pituitary tumor-transforming 1 (PTTG1), transcript variant 2	REVIEWED	mRNA
HCN2	NM_001194.3	7.540	3.60e-41	hyperpolarization activated cyclic nucleotide gated potassium and sodium channel 2 (HCN2)	REVIEWED	mRNA
GDF1	NM_001492.5	4.530	3.60e-41	growth differentiation factor 1 (GDF1)	REVIEWED	mRNA
CERS1	NM_021267.4	4.530	3.60e-41	ceramide synthase 1 (CERS1), transcript variant 1	REVIEWED	mRNA
PLLP	NM_015993.2	8.110	5.92e-41	plasmolipin (PLLP)	VALIDATED	mRNA
TF	NM_001063.3	10.000	6.75e-41	transferrin (TF)	REVIEWED	mRNA
ENDOD1	NM_015036.2	6.100	9.21e-41	endonuclease domain containing 1 (ENDOD1)	VALIDATED	mRNA
HSPA2	NM_021979.3	7.390	4.15e-40	heat shock protein family A (Hsp70) member 2 (HSPA2)	VALIDATED	mRNA
CCNB2	NM_004701.3	-5.730	9.88e-39	cyclin B2 (CCNB2)	REVIEWED	mRNA
PLXNB3	NM_005393.2	7.370	2.34e-38	plexin B3 (PLXNB3), transcript variant 1	REVIEWED	mRNA
ZNF365	NM_014951.2	6.070	7.74e-37	zinc finger protein 365 (ZNF365), transcript variant A	REVIEWED	mRNA
SPNS2	NM_001124758.2	5.330	5.54e-36	sphingolipid transporter 2 (SPNS2)	REVIEWED	mRNA
SNORA109	NR_132964.1	5.890	3.82e-35	small nucleolar RNA, H/ACA box 109 (SNORA109)	VALIDATED	small nucleolar RNA
CNP	NM_033133.4	6.210	1.82e-33	2’,3’-cyclic nucleotide 3’ phosphodiesterase (CNP), transcript variant 1	VALIDATED	mRNA
MT3	NM_005954.2	5.320	3.06e-30	metallothionein 3 (MT3)	VALIDATED	mRNA
HNRNPA1	NM_002136.3	-2.770	4.07e-30	heterogeneous nuclear ribonucleoprotein A1 (HNRNPA1), transcript variant 1	REVIEWED	mRNA
HNRNPA0	NM_006805.3	-2.720	1.04e-29	heterogeneous nuclear ribonucleoprotein A0 (HNRNPA0)	REVIEWED	mRNA
GPR37	NM_005302.3	7.320	2.04e-29	G protein-coupled receptor 37 (GPR37)	REVIEWED	mRNA
QSER1	NM_001076786.2	-3.210	2.75e-29	glutamine and serine rich 1 (QSER1)	VALIDATED	mRNA
KIF1C	XM_005256424.2	5.500	2.89e-29	PREDICTED: Homo sapiens kinesin family member 1C (KIF1C), transcript variant X1	MODEL	mRNA
MAD2L1	NM_002358.3	-4.170	5.12e-29	mitotic arrest deficient 2 like 1 (MAD2L1)	REVIEWED	mRNA
LINC01158	NR_131235.1	-5.810	5.43e-29	long intergenic non-protein coding RNA 1158 (LINC01158), transcript variant 4	VALIDATED	long non-coding RNA
APLP1	NM_005166.4	4.200	2.78e-28	amyloid beta precursor like protein 1 (APLP1), transcript variant 2	VALIDATED	mRNA
TMEM161B-AS1	NR_105020.1	-5.460	8.86e-28	TMEM161B antisense RNA 1 (TMEM161B-AS1), transcript variant 5	VALIDATED	long non-coding RNA
PTGFRN	XM_017001874.1	-3.050	9.48e-28	PREDICTED: Homo sapiens prostaglandin F2 receptor inhibitor (PTGFRN), transcript variant X1	MODEL	mRNA
CKS2	NM_001827.2	-5.200	1.70e-27	CDC28 protein kinase regulatory subunit 2 (CKS2)	REVIEWED	mRNA
TSC22D4	NR_130119.1	6.340	2.23e-27	TSC22 domain family member 4 (TSC22D4), transcript variant 4	VALIDATED	non-coding RNA
TKTL1	NM_001145934.1	-6.410	2.37e-27	transketolase like 1 (TKTL1), transcript variant 3	REVIEWED	mRNA
DCHS1	NM_003737.3	-2.880	4.10e-27	dachsous cadherin-related 1 (DCHS1)	REVIEWED	mRNA
TACC3	NM_006342.2	-4.450	9.17e-27	transforming acidic coiled-coil containing protein 3 (TACC3)	REVIEWED	mRNA
CENPH	NM_022909.3	-4.570	1.34e-26	centromere protein H (CENPH)	REVIEWED	mRNA
ENPP4	NM_014936.4	4.580	3.24e-26	ectonucleotide pyrophosphatase/phosphodiesterase 4 (putative) (ENPP4)	VALIDATED	mRNA
MEX3B	NM_032246.4	-4.320	3.24e-26	mex-3 RNA binding family member B (MEX3B)	REVIEWED	mRNA
PRNP	NM_001080123.2	3.870	3.37e-26	prion protein (PRNP), transcript variant 5	REVIEWED	mRNA

Correlation plot for all genes

This correlation plot uses all the transcripts annotated as PROVISIONAL, REVIEWED or VALIDATED in RefSeq (the high confidence categories)

cor_tlist <- ttg %>% filter(status %in% c('PROVISIONAL', 'REVIEWED', 'VALIDATED')) %>% .[['target_id']]
all <- counts(ddsWald, normalized = TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id')# %>%
up <- all %>%
    filter(target_id %in% (resultsChild %>% as.tibble(.) %>% rownames_to_column(var='target_id') %>% filter(log2FoldChange > 0))[['target_id']])
down <- all %>%
    filter(target_id %in% (resultsChild %>% as.tibble(.) %>% rownames_to_column(var='target_id') %>% filter(log2FoldChange < 0))[['target_id']])
n <- data.frame(All = dim(all)[1], Upregulated = dim(up)[1], Downregulated = dim(down)[1])
kable(n)

All	Upregulated	Downregulated
21848	11382	10466

All transcripts

obs <- 
  counts(ddsWald, normalized = TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>% 
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(est_counts)) %>%
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="grp") %>%
  t() %>%
  as.tibble() #%>%
  #filter_at(vars(matches('Child|Adult')), all_vars(. >= 1))
colnames(obs) <- metadata1 %>% 
  dplyr::select(grp) %>% 
  group_by(grp) %>% 
  summarise(group_count=n()) %>% 
  mutate(group = paste0(grp,'(',group_count,')')) %>%
  .[['group']]
colpal=brewer.pal(n=8, name="RdBu")
corrplot(cor(obs),
         method="circle",
         type="lower",
         col=colpal,
         cl.lim=c(0,1),
         number.cex = 0.6,
         addCoef.col = "black",
         tl.col="black",
         tl.cex = 0.75,
         tl.srt=45,
         diag=FALSE)

Up-regulated

obs <- 
  counts(ddsWald, normalized = TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% (resultsChild %>% as.tibble(.) %>% rownames_to_column(var='target_id') %>% filter(log2FoldChange > 0))[['target_id']]) %>%
  gather(colnames(assay(ddsWald)), key = 'sample', value = 'est_counts') %>% 
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(est_counts)) %>%
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="grp") %>%
  t() %>%
  as.tibble() %>%
  filter_at(vars(matches('Child|Adult')), all_vars(. >= 1))
colnames(obs) <- metadata1 %>% 
  dplyr::select(grp) %>% 
  group_by(grp) %>% 
  summarise(group_count=n()) %>% 
  mutate(group = paste0(grp,'(',group_count,')')) %>%
  .[['group']]
colpal=brewer.pal(n=8, name="RdBu")
corrplot(cor(obs),
         method="circle",
         type="lower",
         col=colpal,
         #cl.lim=c(0,1),
         number.cex = 0.6,
         addCoef.col = "black",
         tl.col="black",
         tl.cex = 0.75,
         tl.srt=45,
         diag=FALSE)

Down-regulated

obs <- 
  counts(ddsWald, normalized=TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% (resultsChild %>% as.tibble(.) %>% rownames_to_column(var='target_id') %>% filter(log2FoldChange < 0))[['target_id']]) %>%
  gather(colnames(assay(ddsWald)), key = 'sample', value = 'est_counts') %>% 
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(est_counts)) %>%
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="grp") %>%
  t() %>%
  as.tibble() %>%
  filter_at(vars(matches('Child|Adult')), all_vars(. >= 1))
colnames(obs) <- metadata1 %>% 
  dplyr::select(grp) %>% 
  group_by(grp) %>% 
  summarise(group_count=n()) %>% 
  mutate(group = paste0(grp,'(',group_count,')')) %>%
  .[['group']]
colpal=brewer.pal(n=8, name="RdBu")
corrplot(cor(obs),
         method="circle",
         type="lower",
         col=colpal,
         #cl.lim=c(0,1),
         number.cex = 0.6,
         addCoef.col = "black",
         tl.col="black",
         tl.cex = 0.75,
         tl.srt=45,
         diag=FALSE)

Correlation plot for “developmental” genes

This correlation plot uses the ‘significant’ transcripts, which are different between the control child samples and the fetal samples.

all <- counts(ddsWald, normalized = TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  # filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% signifChild$target_id)
up <- all %>%
  filter(target_id %in% (signifChild %>% filter(log2FoldChange > 0))[['target_id']])
down <- all %>%
  filter(target_id %in% (signifChild %>% filter(log2FoldChange < 0))[['target_id']])
n <- data.frame(All = dim(all)[1], Upregulated = dim(up)[1], Downregulated = dim(down)[1])
kable(n)

All	Upregulated	Downregulated
1942	1174	768

All transcripts

obs <- 
  counts(ddsWald, normalized = TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% signifChild$target_id) %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>% 
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(est_counts)) %>%
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="grp") %>%
  t() %>%
  as.tibble() %>%
  filter_all(all_vars(. >= 0))
colnames(obs) <- metadata1 %>%
  dplyr::select(grp) %>%
  group_by(grp) %>%
  summarise(group_count=n()) %>%
  mutate(group = paste0(grp,'(',group_count,')')) %>%
  .[['group']]
colpal=brewer.pal(n=8, name="RdBu")
corrplot(cor(obs, method='pearson'),
         method="circle",
         type="lower",
         col=colpal,
         #cl.lim=c(0,1),
         number.cex = 0.6,
         addCoef.col = "black",
         tl.col="black",
         tl.cex = 0.75,
         tl.srt=45,
         diag=FALSE)

Up-regulated

obs <- 
  counts(ddsWald, normalized = TRUE) %>%
  # assay(rld) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% (signifChild %>% filter(log2FoldChange > 0))[['target_id']]) %>%
  gather(colnames(assay(ddsWald)), key = 'sample', value = 'est_counts') %>% 
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = log(median(est_counts)+0.01)) %>%
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="grp") %>%
  t() %>%
  as.tibble() %>%
  filter_all( all_vars(. >= 0))
colnames(obs) <- metadata1 %>% 
  dplyr::select(grp) %>% 
  group_by(grp) %>% 
  summarise(group_count=n()) %>% 
  mutate(group = paste0(grp,'(',group_count,')')) %>%
  .[['group']]
colpal=brewer.pal(n=8, name="RdBu")
corrplot(cor(obs),
         method="circle",
         type="lower",
         col=colpal,
         #cl.lim=c(0,1),
         number.cex = 0.6,
         addCoef.col = "black",
         tl.col="black",
         tl.cex = 0.75,
         tl.srt=45,
         diag=FALSE)

Down-regulated

obs <- 
  counts(ddsWald, normalized=TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% (signifChild %>% filter(log2FoldChange < 0))[['target_id']]) %>%
  gather(colnames(assay(ddsWald)), key = 'sample', value = 'est_counts') %>% 
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(est_counts)) %>%
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="grp") %>%
  t() %>%
  as.tibble() %>%
  filter_at(vars(matches('Child|Adult')), all_vars(. >= 1))
colnames(obs) <- metadata1 %>% 
  dplyr::select(grp) %>% 
  group_by(grp) %>% 
  summarise(group_count=n()) %>% 
  mutate(group = paste0(grp,'(',group_count,')')) %>%
  .[['group']]
colpal=brewer.pal(n=8, name="RdBu")
corrplot(cor(obs),
         method="circle",
         type="lower",
         col=colpal,
         #cl.lim=c(0,1),
         number.cex = 0.6,
         addCoef.col = "black",
         tl.col="black",
         tl.cex = 0.75,
         tl.srt=45,
         diag=FALSE)

Estimated counts distributions

All genes

Unnormalized

counts(ddsWald, normalized=FALSE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  ggplot(aes(x=log2(est_counts), color=sample)) + geom_density() +
  geom_vline(aes(xintercept=0.5), color='red') +
  scale_x_continuous(breaks = round(seq(0, 20, by = 0.5),1)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

Normalized

counts(ddsWald, normalized=TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  ggplot(aes(x=log2(est_counts), color=sample)) + geom_density() +
  geom_vline(aes(xintercept=0.5), color='red') +
  scale_x_continuous(breaks = round(seq(0, 20, by = 0.5),1)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

For “developmental” genes

Unnormalized

counts(ddsWald, normalized=FALSE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% signifChild$target_id) %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  ggplot(aes(x=log2(est_counts), color=sample)) + geom_density() +
  geom_vline(aes(xintercept=0.5), color='red') +
  scale_x_continuous(breaks = round(seq(0, 20, by = 0.5),1)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

Normalized

counts(ddsWald, normalized=TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% signifChild$target_id) %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  ggplot(aes(x=log2(est_counts), color=sample)) + geom_density() +
  geom_vline(aes(xintercept=0.5), color='red') +
  scale_x_continuous(breaks = round(seq(0, 20, by = 0.5),1)) +
  theme(axis.text.x = element_text(angle = 90, hjust = 1))

Correlation plot for “developmental” genes ranked by “closeness” to fetal state

This correlation plot uses the “developmental” genes split in quartiles of closeness to fetal state

all_q <- counts(ddsWald, normalized = TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% signifChild[['target_id']]) %>%
  inner_join((resultsChildArev %>% as_tibble() %>% rownames_to_column(var='target_id')), by='target_id') %>%
  arrange(pvalue) %>% 
  mutate(direction = factor(if_else(target_id %in% (signifChild %>% filter(log2FoldChange > 0))[['target_id']], 'up', 'down'), levels=c('up','down')))
kable(all_q %>% group_by(direction) %>% summarise(count=n()))

direction	count
up	1032
down	525

autism_close_tlist <- (all_q %>% filter(pvalue < 0.99))[['target_id']]

Heatmap of developmental genes

Child

myttg <- ttg %>% 
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% autism_close_tlist) %>%
  dplyr::select(c('target_id', 'symbol'))
symbols <- myttg$symbol
names(symbols) <- myttg$target_id
sel_vars <- metadata1 %>%
  filter(age == 'Child' | grepl('SVZ', sample)) %>%
  dplyr::select('sample') %>%
  as.list() %>% unlist() %>% as.vector()
obs <-  
  counts(ddsWald, normalized=TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% autism_close_tlist) %>%
  dplyr::select(target_id, sel_vars) %>%
  gather(sel_vars, key = 'sample', value = 'est_counts') %>%
  mutate(est_counts=log2(est_counts+1)) %>% 
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="sample") %>%
  as.matrix()
colnames(obs) <- symbols[colnames(obs)]
labels <- metadata1 %>%
  filter(age == 'Child' | grepl('SVZ', sample)) %>%
  dplyr::select(c('sample','age')) %>%
  arrange(sample) %>%
  remove_rownames() %>%
  column_to_rownames(var='sample')
Rowv <- obs %>%
  dist(method="euclidean") %>%
  hclust(method = "complete") %>%
  as.dendrogram %>%
  # rotate(c(6:7,1:5)) %>%
  color_branches(k = 3) %>%
  set("branches_lwd", 4)
obs[which(obs == 0)] <- 1
par(mar=c(5,4,4,2)+0.1)
heatmap.2(obs, Rowv = Rowv, trace="none", margins=c(5,9), srtCol = 45,
          symm = F, symkey = F, symbreaks = F, scale="none")

Adult

myttg <- ttg %>% 
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% autism_close_tlist) %>%
  dplyr::select(c('target_id', 'symbol'))
symbols <- myttg$symbol
names(symbols) <- myttg$target_id
sel_vars <- metadata1 %>%
  filter(age == 'Adult' | grepl('SVZ', sample)) %>%
  dplyr::select('sample') %>%
  as.list() %>% unlist() %>% as.vector()
obs <-
  counts(ddsWald, normalized=TRUE) %>%
  as.data.frame() %>%
  rownames_to_column('target_id') %>%
  filter(target_id %in% cor_tlist) %>%
  filter(target_id %in% autism_close_tlist) %>%
  dplyr::select(target_id, sel_vars) %>%
  gather(sel_vars, key = 'sample', value = 'est_counts') %>%
  mutate(est_counts=log2(est_counts+1)) %>% 
  spread(target_id, est_counts) %>%
  remove_rownames() %>%
  column_to_rownames(var="sample") %>%
  as.matrix()
colnames(obs) <- symbols[colnames(obs)]
labels <- metadata1 %>%
  filter(age == 'Adult' | grepl('SVZ', sample)) %>%
  dplyr::select(c('sample','age')) %>%
  arrange(sample) %>%
  remove_rownames() %>%
  column_to_rownames(var='sample')
Rowv <- obs %>%
  dist(method="euclidean") %>%
  hclust(method = "complete") %>%
  as.dendrogram %>%
  # rotate(c(6:7,1:5)) %>%
  color_branches(k = 3) %>%
  set("branches_lwd", 4)
obs[which(obs == 0)] <- 1
par(mar=c(5,4,4,2)+0.1)
heatmap.2(obs, Rowv = Rowv, trace="none", margins=c(5,9), srtCol = 45,
          symm = F, symkey = F, symbreaks = F, scale="none")

Scatterplots

Child and Adult groups vs Fetal

obs <-
  all_q %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(log2(est_counts+0.001))) %>%
  spread(target_id, est_counts) %>%
  column_to_rownames(var='grp') %>%
  t()
rn <- rownames(obs)
obs <- as.tibble(obs)
obs[['target_id']] <- rn
obs <- 
  obs %>% 
  gather('ChildControl','ChildAutism','AdultControl', 'AdultAutism', key='Group', value='est_counts') %>%
  mutate(Group = factor(Group, levels=c('ChildControl','ChildAutism','AdultControl', 'AdultAutism')))
obs %>%
  ggplot(aes(x=FetalControl)) +
  geom_hex(aes(y=est_counts)) +
  geom_abline(aes(intercept=0, slope=1), color='red') +
  geom_abline(aes(intercept=-1, slope=1), color=muted('red')) +
  geom_abline(aes(intercept=+1, slope=1), color=muted('red')) +
  labs(y = expression(log2('estimated counts'))) +
  facet_wrap(~Group, nrow=2, ncol=2)

ChildAutism vs Fetal by sample

obs1 <-
  all_q %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  full_join(metadata1,  by = 'sample') %>%
  dplyr::select(c('grp','target_id','est_counts')) %>%
  group_by(grp, target_id) %>%
  summarise(est_counts = median(log2(est_counts+0.001))) %>%
  spread(target_id, est_counts) %>%
  column_to_rownames(var='grp') %>%
  t()
rn <- rownames(obs1)
obs1 <- as.tibble(obs1) %>% select('FetalControl')
obs1[['target_id']] <- rn
obs2 <-
  all_q %>%
  gather(colnames(ddsWald), key = 'sample', value = 'est_counts') %>%
  full_join(metadata1,  by = 'sample') %>%
  filter(grp == 'ChildAutism') %>%
  mutate(SampleID=sample) %>%
  dplyr::select(c('SampleID', 'target_id', 'est_counts'))
obs <- 
  obs1 %>% full_join(obs2, by='target_id')
obs %>%
  ggplot(aes(x=FetalControl)) +
  geom_point(aes(y=log2(est_counts+0.01), color=as.factor(SampleID)), alpha=0.6) +
  geom_abline(aes(intercept=0, slope=1), color='red') +
  geom_abline(aes(intercept=-1, slope=1), color=muted('red')) +
  geom_abline(aes(intercept=+1, slope=1), color=muted('red')) +
  labs(y = expression(log2('estimated counts'))) +
  scale_color_brewer(type='qual', palette = 'Dark2', name = "SampleID") + 
  facet_wrap(~SampleID, nrow=2, ncol=2)

Gene trajectories

Up-regulated

for (i in  (all_q %>% filter(pvalue < 0.99 & direction == 'up') )[['target_id']]) {
  obs1 <- txi.salmon$counts %>%
    as.data.frame() %>%
    rownames_to_column('target_id') %>%
    filter(target_id == i) %>%
    gather(colnames(ddsWald), key = 'sample', value = 'est_counts_raw') %>%
    left_join(metadata1, by = 'sample') %>%
    left_join(ttg, by = 'target_id')
  obs2 <- txi.salmon$variance %>%
    as.data.frame() %>%
    rownames_to_column('target_id') %>%
    filter(target_id == i) %>%
    gather(colnames(ddsWald), key = 'sample', value = 'variance_raw')
  obs3 <- normalizationFactors(ddsWald) %>%
    as.data.frame() %>%
    rownames_to_column('target_id') %>%
    filter(target_id == i) %>%
    gather(colnames(ddsWald), key = 'sample', value = 'nf')
  obs1$grp <- factor(obs1$grp, levels = c('FetalControl', 'ChildAutism', 'ChildControl', 'AdultAutism', 'AdultControl'))
  obs <- obs1 %>%
    left_join(obs2) %>%
    left_join(obs3) %>%
    mutate(est_counts = est_counts_raw/nf,
           est_counts_lower = est_counts_raw/nf - 2*sqrt(variance_raw/nf),
           est_counts_upper = est_counts_raw/nf + 2*sqrt(variance_raw/nf)) %>%
    arrange(grp, est_counts)
  obs$sample=factor(obs$sample, levels=obs$sample)
  symbol <- obs$symbol[1]
  cat('####', symbol, '\n\n')
  g <- ggplot(obs, aes(x=sample, y=log(est_counts + 0.01), ymin=log(est_counts_lower + 0.01), ymax=log(est_counts_upper + 0.01), color=grp)) +
    geom_crossbar() +
    labs(title = paste0(obs$symbol[1],' (',obs$target_id[1],')'),
         subtitle = str_c(str_wrap(obs$description[1], width=55), collapse="\n"),
         x = 'sample',
         y = expression(log('estimated counts')),
         fill = "Group") +
    theme(legend.title = element_text(face = "bold.italic"),
          legend.box.background = element_rect(colour = "black", fill = alpha('white', alpha = 0.6)),
          axis.text.x = element_text(angle = 45, hjust = 1),
          plot.title = element_text(hjust = 0.4),
          plot.subtitle = element_text(hjust = 0.4)
    )
  print(g)
  cat('\n\n')
}

P2RY12

ARHGDIB

BTG2

TBX18

IL16

VWF

MPEG1

SOCS3

F5

CYBB

SLC4A2

PCOLCE

TLN1

ZFP36

PROCR

PURA

CPD

FAM129A

C2orf40

NFKBIA

MYH9

INPP4B

GSTM4

TCN2

MEG3

FOS

MCAM

STRN

ADI1

CSF1R

SLC18A2

ARHGAP4

EGR1

SOX17

KLF4

PLXND1

CLN8

PDGFB

C10orf10

SUN2

ESRG

LAPTM5

TPD52L2

Down-regulated

for (i in (all_q %>% filter(pvalue < 0.99 & direction == 'down') )[['target_id']]) {
  obs1 <- txi.salmon$counts %>%
    as.data.frame() %>%
    rownames_to_column('target_id') %>%
    filter(target_id == i) %>%
    gather(colnames(ddsWald), key = 'sample', value = 'est_counts_raw') %>%
    left_join(metadata1, by = 'sample') %>%
    left_join(ttg, by = 'target_id')
  obs2 <- txi.salmon$variance %>%
    as.data.frame() %>%
    rownames_to_column('target_id') %>%
    filter(target_id == i) %>%
    gather(colnames(ddsWald), key = 'sample', value = 'variance_raw')
  obs3 <- normalizationFactors(ddsWald) %>%
    as.data.frame() %>%
    rownames_to_column('target_id') %>%
    filter(target_id == i) %>%
    gather(colnames(ddsWald), key = 'sample', value = 'nf')
  obs1$grp <- factor(obs1$grp, levels = c('FetalControl', 'ChildAutism', 'ChildControl', 'AdultAutism', 'AdultControl'))
  obs <- obs1 %>%
    left_join(obs2) %>%
    left_join(obs3) %>%
    mutate(est_counts = est_counts_raw/nf,
           est_counts_lower = est_counts_raw/nf - 2*sqrt(variance_raw/nf),
           est_counts_upper = est_counts_raw/nf + 2*sqrt(variance_raw/nf)) %>%
    arrange(grp, -est_counts)
  obs$sample=factor(obs$sample, levels=obs$sample)
  symbol <- obs$symbol[1]
  cat('####', symbol, '\n\n')
  g <- ggplot(obs, aes(x=sample, y=log(est_counts + 0.01), ymin=log(est_counts_lower + 0.01), ymax=log(est_counts_upper + 0.01), color=grp)) +
    geom_crossbar() +
    labs(title = paste0(obs$symbol[1],' (',obs$target_id[1],')'),
         subtitle = str_c(str_wrap(obs$description[1], width=55), collapse="\n"),
         x = 'sample',
         y = expression(log('estimated counts')),
         fill = "Group") +
    theme(legend.title = element_text(face = "bold.italic"),
          legend.box.background = element_rect(colour = "black", fill = alpha('white', alpha = 0.6)),
          axis.text.x = element_text(angle = 45, hjust = 1),
          plot.title = element_text(hjust = 0.4),
          plot.subtitle = element_text(hjust = 0.4)
    )
  print(g)
  cat('\n\n')
}

Autism RNASeq correlation analysis

Razvan Sultana

2018-04-13

Preparation

Load packages and set options

Import metadata

SVZ tissue analysis

DESeq2 analysis

Select SVZ metadata

Read RefSeq transcript-to-genes data

Read salmon output files

Perform Wald test

DESeq results

Control - Fetal

Autism - Fetal

QC

PCA

Volcano plots

Control - Fetal

ChildControl - Fetal

AdultControl - Fetal

Autism - Fetal

ChildAutism - Fetal

AdultAutism - Fetal

Gene lists

Control - Fetal

ChildControl - Fetal

AdultControl - Fetal

Autism - Fetal

ChildAutism - Fetal

AdultAutism - Fetal

Correlation plot for all genes

All transcripts

Up-regulated

Down-regulated

Correlation plot for “developmental” genes

All transcripts

Up-regulated

Down-regulated

Estimated counts distributions

All genes

Unnormalized

Normalized

For “developmental” genes

Unnormalized

Normalized

Correlation plot for “developmental” genes ranked by “closeness” to fetal state

Heatmap of developmental genes

Child

Adult

Scatterplots

Child and Adult groups vs Fetal

ChildAutism vs Fetal by sample

Gene trajectories

Up-regulated

P2RY12

ARHGDIB

BTG2

TBX18

IL16

VWF

MPEG1

SOCS3

F5

CYBB

SLC4A2

PCOLCE

TLN1

ZFP36

PROCR

PURA

CPD

FAM129A

C2orf40

NFKBIA

MYH9

INPP4B

GSTM4

TCN2

MEG3