{"type":"doc","content":[{"type":"paragraph","content":[{"text":"Legacy","type":"text"}]},{"type":"codeBlock","content":[{"text":"##STEP1: import gene counts\ncountData <- read.table('salmon.merge.gene_counts.tsv', header = TRUE, sep = \"\\t\")\nhead(countData)\n\n##STEP2: import metaData\nmetaData <- read.csv('metadata.tsv*', header = TRUE, sep = \",\")\nmetaData\n\n##STEP3: Construct DESEQDataSet Object\ndds <- DESeqDataSetFromMatrix(countData=countData, \n colData=metaData, \n design=~celltype, tidy = TRUE)\n\n#define the reference group set as the first mentioned in levels\ndds$celltype <- factor(dds$celltype, levels = c(\"CON\",\"TB\"))\n\n##STEP4: Run DESEQ function\ndds <- DESeq(dds)\n\n\n##STPE5: Take a look at the results table\nres <- results(dds)\nhead(results(dds, tidy=TRUE)) #let's look at the results table\n\n\n##STEP6: Summary of differential gene expression\nsummary(res) #summary of results\n\n\n##STEP7: Sort summary list by p-value\nres <- res[order(res$padj),]\nhead(res)\nhead(res, 9)\n\n##STEP8: save DEGs to file\nwrite.table(res, file = \"NoHeadInjury_CON__vs__HeadInjury_TB__DESeq2_toptags_DEGs.txt\")\n\n\n##STEP9: save normalised counts to file\nfoo <- counts(dds, normalized = TRUE)\nwrite.csv(foo, file=\"NoHeadInjury_CON__vs__HeadInjury_TB__DESeq2_toptags_DEGs__Normalised_Counts.csv\")\n\n\n##STEP10: PCA\n#version 1\npdf(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__PCA.pdf\")\nvsdata <- vst(dds, blind=TRUE, nsub = 100)\nplotPCA(vsdata, intgroup=\"dex\") \ndev.off()\n\n#version 2:\npdf(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB_label.pdf\")\np <- plotPCA(vsdata, intgroup=\"dex\") \np + geom_text_repel(aes(label=substr(name, start = 1, stop = 5)),size = 2, min.segment.length = 0, seed = 42, box.padding = 0.5, max.overlaps = 20)\ndev.off()\n\n## jpeq\nlibrary(Cairo)\n\nCairo(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__PCA_label.jpeg\",\n type=\"png\",\n units=\"px\", \n width=2500, \n height=2000, \n pointsize=12, \n dpi=300)\np <- plotPCA(vsdata, intgroup=\"dex\") \np + geom_text_repel(aes(label=substr(name, start = 1, stop = 5)),size = 2, min.segment.length = 0, seed = 42, box.padding = 0.5, max.overlaps = 20)\ndev.off()\n\n\n#STEP11: VOLCANO PLOTS\n#basic volcano plot \nEnhancedVolcano(res,\n lab = rownames(res),\n x = 'log2FoldChange',\n y = 'pvalue')\n\n#save enhanced volcano plot to file\npdf(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__enhanced_volcano_plot_2.0FC_p1e-2.pdf\")\nEnhancedVolcano(res,\n lab = rownames(res),\n x = 'log2FoldChange',\n y = 'pvalue',\n title = 'NoHeadInjury_CON__vs__HeadInjury_TB',\n pCutoff = 1e-2,\n FCcutoff = 1.5,\n pointSize = 1.0,\n labSize = 3.5)\ndev.off()\n\npdf(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__enhanced_volcano_plot_2.0FC_p10e-3.pdf\")\nEnhancedVolcano(res,\n lab = rownames(res),\n x = \"log2FoldChange\",\n y = \"pvalue\",\n pCutoff = 10e-3,\n FCcutoff = 2,\n ylim = c(0, -log10(10e-12)),\n pointSize = c(ifelse(res$log2FoldChange>2, 2, 2)),\n labSize = 3.0,\n #shape = c(6, 6, 19, 16),\n title = \"DESeq2: NoHeadInjury_CON vs. HeadInjury_TB\",\n subtitle = \"Differential expression\",\n caption = bquote(~Log[2]~ \"fold change cutoff, 2; p-value cutoff, 10e-3\"),\n #legendPosition = \"right\",\n legendLabSize = 8,\n col = c(\"grey30\", \"forestgreen\", \"royalblue\", \"red2\"),\n colAlpha = 0.9,\n drawConnectors = TRUE,\n hline = c(10e-3),\n widthConnectors = 0.5)\ndev.off()\n\n\n#save as high resolution jpeg\nlibrary(Cairo)\n\n## Initialize the device for saving first plot.\n## Plot will have length=500 pixel, width=500 pixel\n\nCairo(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__enhanced_volcano_plot_2.0FC_p1e-3.jpeg\",\n type=\"png\",\n units=\"px\", \n width=2800, \n height=2500, \n pointsize=12, \n dpi=300)\nEnhancedVolcano(res,\n lab = rownames(res),\n x = 'log2FoldChange',\n y = 'pvalue',\n title = 'NoHeadInjury_CON__vs__HeadInjury_TB',\n pCutoff = 1e-3,\n FCcutoff = 2.0,\n pointSize = 1.0,\n labSize = 3.0)\ndev.off()\n\n\n#option2:\nCairo(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__enhanced_volcano_plot_2.0FC_p10e-3_v2.jpeg\",\n type=\"png\",\n units=\"px\", \n width=2500, \n height=2200, \n pointsize=14, \n dpi=300)\nEnhancedVolcano(res,\n lab = rownames(res),\n x = \"log2FoldChange\",\n y = \"pvalue\",\n pCutoff = 10e-3,\n FCcutoff = 2,\n ylim = c(0, -log10(10e-12)),\n pointSize = c(ifelse(res$log2FoldChange>2, 2, 2)),\n labSize = 3.5,\n #shape = c(6, 6, 19, 16),\n title = \"DESeq2: NoHeadInjury_CON vs. HeadInjury_TB\",\n subtitle = \"Differential expression\",\n caption = bquote(~Log[2]~ \"fold change cutoff, 2; p-value cutoff, 10e-3\"),\n #legendPosition = \"right\",\n legendLabSize = 8,\n col = c(\"grey30\", \"forestgreen\", \"royalblue\", \"red2\"),\n colAlpha = 0.9,\n drawConnectors = TRUE,\n hline = c(10e-3),\n widthConnectors = 0.5)\ndev.off()\n\n\n\n\n#save as high res tiff\n\ntiff(file=\"Figure_NoHeadInjury_CON__vs__HeadInjury_TB__enhanced_volcano_plot_2.0FC_p1e-5.tiff\",\n width=8, height=8, units=\"in\", res=150)\nEnhancedVolcano(res,\n lab = rownames(res),\n x = 'log2FoldChange',\n y = 'pvalue',\n title = 'NoHeadInjury_CON__vs__HeadInjury_TB',\n pCutoff = 1e-5,\n FCcutoff = 2.0,\n pointSize = 1.5,\n labSize = 3.5)\n\ndev.off()","type":"text"}]}],"version":1}

Browser not supported