Question: Identifying variants differing between control/treatment
0
gravatar for rmf
23 months ago by
rmf760
rmf760 wrote:

I have a single multisample (48 samples) sample VCF file generated by GATK joint genotyping. I have further filtered by site level quality and such. Of these 48 samples, 24 samples are control and 24 samples are treatment. I am not so much interested between my samples vs reference. I am interested in variants that differ between my control and treated samples. How do I go about doing this? What is the workflow? What tools should I use?

snp indel gatk variant-calling • 671 views
ADD COMMENTlink modified 23 months ago by Pierre Lindenbaum124k • written 23 months ago by rmf760
1
gravatar for Pierre Lindenbaum
23 months ago by
France/Nantes/Institut du Thorax - INSERM UMR1087
Pierre Lindenbaum124k wrote:

using vcfilterjdk http://lindenb.github.io/jvarkit/VcfFilterJdk.html

the following code declare two samples 'cases' S1 and S2, and two samples 'controls' S3 and S5. it then check that all cases have a different genotype than the controls.

 java -jar dist/vcffilterjdk.jar -e 'String cases[]=new String[]{"S1","S2"};String controls[]=new String[]{"S3","S4"}; for(final String S1:cases) {final Genotype G1=variant.getGenotype(S1); if(G1==null ||!G1.isCalled()) continue;for(final String S2:controls) {final Genotype G2=variant.getGenotype(S2); if(G2==null ||!G2.isCalled()) continue;   if(G1.sameGenotype(G2)) return false;}} return true;' input.vcf
ADD COMMENTlink modified 23 months ago • written 23 months ago by Pierre Lindenbaum124k

So, in my case, I just need to replace S1, S2 etc with my sample aliases? Can I source them as a list of sample aliases in a text file?

control.txt

c1
c2
c3

treat.txt

t1
t2
t3
ADD REPLYlink modified 23 months ago • written 23 months ago by rmf760
sure:

$ java -jar dist/vcffilterjdk.jar --body -e 'List<String> cases = null,controls=null;  public Object apply(final VariantContext variant) { if(cases==null) try {cases=IOUtil.slurpLines(new java.io.File("treat.txt")) ; controls= IOUtil.slurpLines(new java.io.File("control.txt")); } catch(Exception e) {throw new RuntimeIOException(e);}; for(final String S1:cases) {final Genotype G1=variant.getGenotype(S1); if(G1==null ||!G1.isCalled()) continue;for(final String S2:controls) {final Genotype G2=variant.getGenotype(S2); if(G2==null ||!G2.isCalled()) continue;   if(G1.sameGenotype(G2)) return false;}} return true;}' input.vcf
ADD REPLYlink written 23 months ago by Pierre Lindenbaum124k

Downloaded, compiled and ran this: java -jar vcffilterjdk.jar --help

Error: A JNI error has occurred, please check your installation and try again
Exception in thread "main" java.lang.NoClassDefFoundError: com/beust/jcommander/ParameterException
        at java.lang.Class.getDeclaredMethods0(Native Method)
        at java.lang.Class.privateGetDeclaredMethods(Class.java:2701)
        at java.lang.Class.privateGetMethodRecursive(Class.java:3048)
        at java.lang.Class.getMethod0(Class.java:3018)
        at java.lang.Class.getMethod(Class.java:1784)
        at sun.launcher.LauncherHelper.validateMainClass(LauncherHelper.java:544)
        at sun.launcher.LauncherHelper.checkAndLoadMain(LauncherHelper.java:526)
Caused by: java.lang.ClassNotFoundException: com.beust.jcommander.ParameterException
        at java.net.URLClassLoader.findClass(URLClassLoader.java:381)
        at java.lang.ClassLoader.loadClass(ClassLoader.java:424)
        at sun.misc.Launcher$AppClassLoader.loadClass(Launcher.java:331)
        at java.lang.ClassLoader.loadClass(ClassLoader.java:357)

I have checked if I have JAVA in my path echo ${PATH} /sw/comp/java/x86_64/sun_jdk1.8.0_92/bin

ADD REPLYlink modified 23 months ago • written 23 months ago by rmf760

hum what is the output of

java -version

please

ADD REPLYlink written 23 months ago by Pierre Lindenbaum124k

Checked all dependencies.

[rf@m1 dist]$ java -version
java version "1.8.0_92"
Java(TM) SE Runtime Environment (build 1.8.0_92-b14)
Java HotSpot(TM) 64-Bit Server VM (build 25.92-b14, mixed mode)

[rf@m1 dist]$ xsltproc --version
Using libxml 20706, libxslt 10126 and libexslt 815
xsltproc was compiled against libxml 20706, libxslt 10126 and libexslt 815
libxslt 10126 was compiled against libxml 20706
libexslt 815 was compiled against libxml 20706

[rf@m1 dist]$ make --version
GNU Make 3.81
Copyright (C) 2006  Free Software Foundation, Inc.

[rf@m1 dist]$ git --version
git version 1.7.1

[rf@m1 dist]$ wget --version
GNU Wget 1.12 built on linux-gnu.
ADD REPLYlink modified 23 months ago • written 23 months ago by rmf760

this is really strange... did you move the jar file ?

ADD REPLYlink written 23 months ago by Pierre Lindenbaum124k

nope.

ADD REPLYlink written 23 months ago by rmf760
1

Deleted all and reinstalled again. It works this time.

ADD REPLYlink written 23 months ago by rmf760
Please log in to add an answer.

Help
Access

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.3.0
Traffic: 1162 users visited in the last hour