bmansfeld / QTLseqr

QTLseqr is an R package for QTL mapping using NGS Bulk Segregant Analysis
64 stars 42 forks source link

error: fewer than one row in the data in runGprimeAnalysis #17

Closed cjchen5 closed 5 years ago

cjchen5 commented 5 years ago

Hi Ben, I got a problem in runGprimeAnalysis when I use version 0.7.5.2 Here is my code, same with your guide.

> df_filt <- runGprimeAnalysis(
+   SNPset = df_filt,
+   windowSize = 1e6,
+   outlierFilter = "deltaSNP")

Here is output:

Counting SNPs in each window...
Calculating tricube smoothed delta SNP index...
Calculating G and G' statistics...
Error in locfit::locfit(Stat ~ locfit::lp(POS, h = windowSize, deg = 0),  : 
  fewer than one row in the data
In addition: There were 37 warnings (use warnings() to see them)
> warnings()
Warning messages:
1: In lfproc(x, y, weights = weights, cens = cens, base = base,  ... :
  Estimated rdf < 1.0; not estimating variance

Here is code I use to get table, same with your guide for -F and -GF:

gatk VariantsToTable -R ${ref}/$reference_file -V ${reference}_Mr_cohort_snps.vcf -F CHROM -F POS -F REF -F ALT -GF AD -GF DP -GF GQ -GF PL -O ${reference}_Mr_db.snps.vcf.table

Here is head of df:

> head(df)
  CHROM   POS REF ALT AD_REF.LOW AD_ALT.LOW DP.LOW GQ.LOW     PL.LOW SNPindex.LOW AD_REF.HIGH AD_ALT.HIGH DP.HIGH GQ.HIGH     PL.HIGH SNPindex.HIGH
1  Chr1 31071   A   G         34         36     70     99  897,0,855    0.5142857          26          22      48      99   522,0,698     0.4583333
2  Chr1 31478   C   T         34         52     86     99 1363,0,844    0.6046512          40          34      74      99  848,0,1099     0.4594595
3  Chr1 33667   A   G         20         48     68     99 1331,0,438    0.7058824          24          29      53      99   765,0,599     0.5471698
4  Chr1 34057   C   T         38         40     78     99 1059,0,996    0.5128205          29          26      55      99   673,0,772     0.4727273
5  Chr1 35239   A   C         25         36     61     99  987,0,630    0.5901639          40          60     100      99 1632,0,1015     0.6000000
6  Chr1 38389   T   C         36         42     78     99 1066,0,906    0.5384615          42          40      82      99  984,0,1105     0.4878049
    REF_FRQ     deltaSNP
1 0.5084746 -0.055952381
2 0.4625000 -0.145191703
3 0.3636364 -0.158712542
4 0.5037594 -0.040093240
5 0.4037267  0.009836066
6 0.4875000 -0.050656660

Could you tell me what might cause that error? Thanks!

bmansfeld commented 5 years ago

Hi, You are showing me the head(df) but running the command on df_filt. df_filt <- runGprimeAnalysis(SNPset = df_filt, windowSize = 1e6, outlierFilter = "deltaSNP") Did you filter snps already? How many rows in the df_filt data frame? Just based on the error msg it seems that the df_filt is empty? Or maybe one of the Chrms is empty. Something is strange with the data. Ben

cjchen5 commented 5 years ago

Hi Ben, Thank you! I guess one of the Chrms is empty after filtering. I will check it later. BTW, I have another question for Gprime plot. The deltaSNP plot from same data is continuous but Gprime plot have a weird gap near the peak. Could you tell me what might cause the gap in Gprime plot?

bmansfeld commented 5 years ago

Hi, The default y axis limits are -0.6 to 0.6, you have a nice peak there going past that range so it's getting cut off. You are probably also getting a warning message about values not being mapped. To remedy this you should add + ylim(-1,1) or some other range after your plotQTLStats() function. Ben

cjchen5 commented 5 years ago

Thanks! But when I change my plot code to: plotQTLStats(SNPset = df_filt, var = "deltaSNP", plotIntervals = TRUE, + ylim(-1,1)) Error in ylim(-1, 1) : could not find function "ylim" Could you please tell me more detail of add + ylim(-1,1)? Thank you very much!

bmansfeld commented 5 years ago

The function is a wrapper for ggplot therefore you need to add + ylim() after the function. I.e. plotQTLStats(....) + ylim(..) Look at the ggplot package to understand more. Ben

cjchen5 commented 5 years ago

Thank you so much! When I library("ggplot2") first and then plotQTLStats(....) + ylim(..) it works.