Markdown: FindMutantVsWildtypeDEGenes-DESeq2.Rmd

Bug:
There was a bug in the MvW markdown. This bug included not containing the accurate results that were seen before in the MvW.34.75-SL5.txt results file. In the MvW-SL5.txt you could see that the statistics were different and had opposite symbols for negative and positive numbers.

Find Bug:

Under 'Perform the Comparisons' section of the markdown the following code was the issue:

rs = results(dds,alpha=Q)

The reason this was the issue was because there was no definition of 'lfcThreshold' so the default was set to 0. Then there was no contrast so the reference wild type, VF36, was being subtracted first and causing incorrect symbols.

Fix Bug:

Now we know what the issue is we can fix the code to show the correct results:

    rs = results(dds,alpha=Q, lfcThreshold = lfcThreshold, 
                 contrast = c("condition", group_2, group_1))

The 'lfcThreshold' is now fixed to 1 which was defined at the beginning of the markdown. So now the correct pvalues should be showing. The contrast now has the 'ARE' (group_2) genotype subtracting with the wiltype 'VF36' (group_1). So now the symbols should be correct.

Note: Q has also been changed back to 0.05 which is what it was originally. Also the comps table, which includes the sum of significant genes, has also changed.

Branch: https://bitbucket.org/mdavis4290/molly5-flavonoid-rnaseq/branch/IGBF-3504
Commit: https://bitbucket.org/mdavis4290/molly5-flavonoid-rnaseq/commits/ef235862778dee176f3e274b1f9eaf6822a2e198

Thanks [~molly]!

Please submit PR when ready.

PR: https://bitbucket.org/hotpollen/flavonoid-rnaseq/pull-requests/37

PR is merged.

Checked new output by looking at MvW-SL5.txt, like so:

cut -f1-9  MvW-SL5.txt | grep V.34.75 | grep A.34.75 | grep Solyc07G000215

The new data file shows test statistic of -16.8, which matches what was obtained from version bba3f03, which is also shown in the first tab of the attached spreadsheet.

Also, the new data file shows the proper sign for the fold-change: negative means the are genotype was less.

Note: we may need to fix any explanatory text in the knitted Markdown about what group1 and group 2 mean. Here, group2 is the numerator and group1 is the denominator in the log2 fold-change.

Testing:

local aloraine$ cut -f1-9  MvW-SL5.txt | grep Solyc07G000215 
Solyc07G000215	A.28	V.28	296	-3.79	0.201	3.39e-77	2.04e-79	-18.9
Solyc07G000215	V.28.15	A.28.15	274	-3.34	0.203	3.96e-27	7.96e-31	-11.5
Solyc07G000215	V.28.30	A.28.30	279	-3.46	0.44	8.96e-06	2.41e-08	-5.58
Solyc07G000215	V.28.45	A.28.45	331	-4.46	0.238	6.27e-44	5.41e-48	-14.6
Solyc07G000215	V.28.75	A.28.75	312	-3.99	0.231	3.03e-34	2.92e-38	-12.9
Solyc07G000215	V.34.15	A.34.15	273	-3.5	0.2	7.28e-32	5.93e-36	-12.5
Solyc07G000215	V.34.30	A.34.30	285	-3.36	0.514	0.00107	4.41e-06	-4.59
Solyc07G000215	V.34.45	A.34.45	323	-4.09	0.212	3.17e-44	3.16e-48	-14.6
Solyc07G000215	V.34.75	A.34.75	339	-4.72	0.222	3.03e-59	3.2e-63	-16.8

Attached barchart showing expression levels. Expression in are genotype is less than VF36 in every pairing or grouping.

Changes needed:

• group1 and group2 need to be fixed for A.28 vs V.28 comparison.

Changed code so that group1 is the numerator and group2 is the denominator.
Also, add lfcThreshold option to A.28 versus V.28 comparison.
Checked with:

local aloraine$ cut -f1-9  MvW-SL5.txt | grep Solyc07G000215 
Solyc07G000215	A.28	V.28	296	-3.79	0.201	1.27e-40	7.41e-44	-13.9
Solyc07G000215	A.28.15	V.28.15	274	-3.34	0.203	3.97e-27	7.99e-31	-11.5
Solyc07G000215	A.28.30	V.28.30	279	-3.46	0.44	8.96e-06	2.41e-08	-5.58
Solyc07G000215	A.28.45	V.28.45	331	-4.46	0.238	6.27e-44	5.41e-48	-14.6
Solyc07G000215	A.28.75	V.28.75	312	-3.99	0.231	3.03e-34	2.92e-38	-12.9
Solyc07G000215	A.34.15	V.34.15	273	-3.5	0.2	7.28e-32	5.93e-36	-12.5
Solyc07G000215	A.34.30	V.34.30	285	-3.36	0.514	0.00107	4.41e-06	-4.59
Solyc07G000215	A.34.45	V.34.45	323	-4.09	0.212	3.17e-44	3.16e-48	-14.6
Solyc07G000215	A.34.75	V.34.75	339	-4.72	0.222	3.03e-59	3.19e-63	-16.8

Notified Muday lab about updated files.
Moving to DONE.