Multiple QTL model analysis with covariate

482 views
Skip to first unread message

jspv2 v

unread,
Jul 29, 2019, 10:01:45 AM7/29/19
to R/qtl discussion
Dear Karl,

Currently, I am performing a QTL analysis with a RIL population (F6)  in R/QTL.
I have identified a covariate which affects the phenotype I use for the mapping, and and using the scanone analysis I have determined that this is most likely an additive covariate.

Now I want to use this information in the Multiple QTL model I am building with the fitqtl and refineqtl options. I was wondering: Do I then add these covariates to the model as follows y~Q1+Q2+Q1*Q1+Covariate,  or do I need to put the model as Q1+Q2+Q1*Q2+  Q1 + Covariate + Q2 + Covariate ?

Also I was thinking,  is it wise to again test here whether I have an additive or interacting covariate? (e.g. adding the covariate as an additive and as interacting covariate to the multiple QTL model, and then compare which model has the best LOD-scores?)

Thanks!

Best Wishes,
Jasper

Karl Broman

unread,
Jul 29, 2019, 1:05:34 PM7/29/19
to rqtl...@googlegroups.com
In the fitqtl/refineqtl functions, you need to create a data frame of
covariates, and you refer to the individual columns by their column names.

If an additive covariate, you'd just write like

  y ~ Q1+Q2+Q1*Q2+Covariate

If you want to include an interaction with QTL 1, you'd do

  y ~ Q1+Q2+Q1*Q2+Covariate + Covariate:Q1

: is for a specific interaction
* includes the interaction plus all marginal effects

Yes I would reevaluate the evidence for the QTL and interactions in the
context of the multiple-QTL model.

karl
> --
> You received this message because you are subscribed to the Google
> Groups "R/qtl discussion" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to rqtl-disc+...@googlegroups.com
> <mailto:rqtl-disc+...@googlegroups.com>.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/rqtl-disc/09239b15-f7e0-4243-a6c6-b4e138eb5a8f%40googlegroups.com
> <https://groups.google.com/d/msgid/rqtl-disc/09239b15-f7e0-4243-a6c6-b4e138eb5a8f%40googlegroups.com?utm_medium=email&utm_source=footer>.

jspv2 v

unread,
Jul 29, 2019, 2:29:20 PM7/29/19
to R/qtl discussion
Dear Karl,

Thanks for your advice!

In both the interval mapping and now in the multiple QTL model i see no evidence for an interacting covariate, so thats nicely in agreement.

One other analysis I now did was scantwo, and without the covariate I find 3 pairs of QTLs. Adding the covariate makes 2 of these pairs not significant anymore in the scantwo analysis.
However, when I try to combine all 3 pairs of QTLs in the multiple QTL model with the covariate, all QTLs have a significant LOD score. I determined the significance treshold using the permutation test of scanone with the covariate added. Is this the right approach? And if so, should I then trust more the results of scantwo, or the multiple QTL model results?

Thanks!

Best Wishes,
Jasper



Op maandag 29 juli 2019 19:05:34 UTC+2 schreef Karl Broman:

Karl Broman

unread,
Jul 29, 2019, 4:43:13 PM7/29/19
to rqtl...@googlegroups.com
I don't fully understand. If you're trying to establish the significance
of QTL pairs/interactions, you need to look at scantwo permutations, not
scanone permutations.

If you're not looking for QTL x covariate interactions, but just want to
include the covariate as an additive term in all models, you'd include
it as an additive covariate in each of scanone and scantwo, and the
permutations for them.

In general, I'd trust the multiple-QTL model over the results of scantwo.

karl
> > an email to rqtl...@googlegroups.com <javascript:>
> > <mailto:rqtl...@googlegroups.com <javascript:>>.
> > To view this discussion on the web visit
> >
> https://groups.google.com/d/msgid/rqtl-disc/09239b15-f7e0-4243-a6c6-b4e138eb5a8f%40googlegroups.com
> <https://groups.google.com/d/msgid/rqtl-disc/09239b15-f7e0-4243-a6c6-b4e138eb5a8f%40googlegroups.com>
>
> >
> <https://groups.google.com/d/msgid/rqtl-disc/09239b15-f7e0-4243-a6c6-b4e138eb5a8f%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/rqtl-disc/09239b15-f7e0-4243-a6c6-b4e138eb5a8f%40googlegroups.com?utm_medium=email&utm_source=footer>>.
>
>
> --
> You received this message because you are subscribed to the Google
> Groups "R/qtl discussion" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to rqtl-disc+...@googlegroups.com
> <mailto:rqtl-disc+...@googlegroups.com>.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/rqtl-disc/d6d77e49-ebff-47ba-a98b-2671e9e52d45%40googlegroups.com
> <https://groups.google.com/d/msgid/rqtl-disc/d6d77e49-ebff-47ba-a98b-2671e9e52d45%40googlegroups.com?utm_medium=email&utm_source=footer>.

jspv2 v

unread,
Jul 29, 2019, 4:54:26 PM7/29/19
to R/qtl discussion
Hi Karl,

Indeed, thats what I did.

So I did a scantwo analysis with and without covariates, and then also determined the significance with the permutation test of scantwo with and without these covariates.
There I find for the scantwo without covariates 3 pairs of QTLs, and for the scantwo.

Then, I tried to build a multiple QTL model (using fitqtl, refineqtl, etc) with all these 3 pairs of QTLs I found with scantwo. In this multiple QTL model, I also added the covariate. Now, I see that the pairs of QTLs that were not significant in the scantwo analysis, appear significant in the multiple QTL model (all LOD score over 3).   For this multiple qtl model, I estimated the significance tresholds using the scanone function with the covariate.
I was wondering if this is the correct way to estimate the significance tresholds for the multiple QTL model?

Op maandag 29 juli 2019 22:43:13 UTC+2 schreef Karl Broman:

Karl Broman

unread,
Jul 29, 2019, 5:07:05 PM7/29/19
to rqtl...@googlegroups.com
My preference is to use penalized LOD scores as described in Manichaikul
et al (2009) https://doi.org/10.1534/genetics.108.094565
This would use scantwo permutations, which I think would be important as
you're three pairs of QTL were identified in scantwo.(vs scanone, as you
write in "For this multiple qtl model, I estimated the significance
tresholds using the scanone function with the covariate.")

karl
> https://groups.google.com/d/msgid/rqtl-disc/d6d77e49-ebff-47ba-a98b-2671e9e52d45%40googlegroups.com
> <https://groups.google.com/d/msgid/rqtl-disc/d6d77e49-ebff-47ba-a98b-2671e9e52d45%40googlegroups.com>
>
> >
> <https://groups.google.com/d/msgid/rqtl-disc/d6d77e49-ebff-47ba-a98b-2671e9e52d45%40googlegroups.com?utm_medium=email&utm_source=footer
> <https://groups.google.com/d/msgid/rqtl-disc/d6d77e49-ebff-47ba-a98b-2671e9e52d45%40googlegroups.com?utm_medium=email&utm_source=footer>>.
>
>
> --
> You received this message because you are subscribed to the Google
> Groups "R/qtl discussion" group.
> To unsubscribe from this group and stop receiving emails from it, send
> an email to rqtl-disc+...@googlegroups.com
> <mailto:rqtl-disc+...@googlegroups.com>.
> To view this discussion on the web visit
> https://groups.google.com/d/msgid/rqtl-disc/288f2988-1766-4cbe-8e55-b93e200ca015%40googlegroups.com
> <https://groups.google.com/d/msgid/rqtl-disc/288f2988-1766-4cbe-8e55-b93e200ca015%40googlegroups.com?utm_medium=email&utm_source=footer>.

jspv2 v

unread,
Aug 2, 2019, 8:26:18 AM8/2/19
to R/qtl discussion
Dear Karl,

Thanks,
I am now using the stepwiseqtl function, with the penalised low-score as described in chapter 9.
Interestingly, with stepwiseQTL, one QTL is detected that I did not detect using either the scanone, scantwo, addqtl, addint functions.
I was wondering what could cause this?  
I do also get several warnings: "In min(diff(a)) : no non-missing arguments to min; returning Inf"
could these warnings mean there is something going wrong with the analysis?

Or are the differences in the significant QTLs more caused by the different ways to calculate the LOD-tresholds? (as for scanone/scantwo the method "em" is used, while for stepwise the penalised LOD score is calculated with the method "imp"?)


Best Wishes,
Jasper


Op maandag 29 juli 2019 23:07:05 UTC+2 schreef Karl Broman:

Karl Broman

unread,
Aug 3, 2019, 10:03:51 AM8/3/19
to rqtl...@googlegroups.com
There are a number of ways that a QTL might be observed from stepwiseqtl that wasn't seen in scanone.
But the em/hk/imp method difference is one possibility. It would be helpful if you showed some more results.
Putative QTL aren't just there or not there, but there's varying degrees of evidence for them.

I can't tell what the min(diff(a)) warnings might be from without further details.

karl

To unsubscribe from this group and stop receiving emails from it, send an email to rqtl-disc+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/rqtl-disc/c437898e-2c06-4553-a0cd-12c4ca004a21%40googlegroups.com.

jspv2 v

unread,
Sep 16, 2019, 3:52:39 PM9/16/19
to R/qtl discussion
Hi Karl,

Thanks!

The min(diff(a)) warnings I get show up after the StepwiseQTL analysis. One other thing I noticed: the pLod values in the StepwiseQTL are negative. What does that mean?
Does it simply mean the QTLs are not significant?   At the scanone analysis, I get only 1 peak that is slightly over the LOD-treshold, the rest is not significant. In the scantwo, I find 1 pair of (different) QTLs that are also just a little bit higher than the significance treshold..

At the stepwiseQTL analysis, I get in total 32 warnings, which all are the same "min(diff(a)" error..

What I also noticed, with the same dataset, I get also twice the "min(diff(a)" error when I use the refineqtl function. Could these errors be related?

Best Wishes,
Jasper





Op zaterdag 3 augustus 2019 16:03:51 UTC+2 schreef Karl Broman:
To unsubscribe from this group and stop receiving emails from it, send an email to rqtl...@googlegroups.com.

Karl Broman

unread,
Sep 16, 2019, 4:02:48 PM9/16/19
to rqtl...@googlegroups.com
A negative pLOD means the model is worse than the null model (with no QTL) which has pLOD = 0.

To diagnose the min(diff(a)) I would need to see some code and ideally example data. I need to be able to reproduce the problem before I can diagnose and fix it.

karl
To unsubscribe from this group and stop receiving emails from it, send an email to rqtl-disc+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/rqtl-disc/e762a708-407a-4e8b-8180-c5addd139b4c%40googlegroups.com.
Message has been deleted

jspv2 v

unread,
Sep 17, 2019, 12:18:05 PM9/17/19
to R/qtl discussion
Hi Karl,

As in my PhD, I'm working together with several companies, the data is confidential, so not sure if I can share the real data here.
But I can share with you the code/script I used:


So I start with the read.cross function, and then the jittermap function:
stuntRIL <- read.cross("csv", file="newmapXX.csv", genotypes = c ("AA", "AB", "BB"), na.strings="-", alleles=c("A", "B"), crosstype="riself")
stuntRIL2 <- jittermap(stuntRIL)

Then I add the covariate I work with:
covarAB <- cbind(Area, Block)
(I use a blocking design in my experiment, so I use the Block as a covariate, and then also size of the individuals as another covariate)

Then I perform the scanone, with and without the covariate:
out.0 <-scanone(stuntRIL2, method="em", pheno.col=1)
perm.0 <- scanone(stuntRIL2, method="em", pheno.col=1, n.perm=500, n.cluster=8)
summary(perm.0, a=c(0.05,0.01))
plot(out.0, ylab="LOD score", ylim=c(0,5))  + abline(h=2.7)

out.AB <-scanone(stuntRIL2, method="em", pheno.col=1, addcovar=covarAB)
perm.AB <- scanone(stuntRIL2, method="em", pheno.col=1, addcovar=covarAB, n.perm=1000, n.cluster=8)
summary(perm.AB, a=c(0.05,0.01))
plot(out.AB, ylab="LOD score", ylim=c(0,5))  + abline(h=2.85)

Until so far, I get no errors. In both with and withouth the covariate, I get several peaks, but only one reaches over the lod-treshold with a top value of around 2.9. So really not much strong evidence for QTLs. That could be because the phenotype I work with has a high degree of variability.

Then I performed the scantwo function:
out2.0 <- scantwo(stuntRIL2, verbose=TRUE, pheno.col=1)
max(out2.0)
summary(out2.0, perms=Perm0, alphas=c(0.05, 0.05, 0.05, 0.05, 0.05), pvalues=TRUE, clean.output=TRUE)

Also there no problems. As results I get one pair of QTLs when I set the permutation treshold at a=0.05 (with LODs just slightly higher than the treshold). I get 2 pairs when I set the threshold at a=0.10.



The problems then start when I start to build the multiple QTL model with the following script:

stuntRIL2 <- sim.geno(RIL2, step=2, n.draws=400, err=0.01)
qtl2.a <- makeqtl(RIL2, chr=c(), pos=c())    (Chr. numbers and positions I cannot publicly post on google due to the confidentiality of my project)
out.fqB2.a <- fitqtl(RIL2, qtl=qtl2.a, formula=y~Q1+Q2+Q3+Q4+Q5+Q6 +Q1*Q2 + Q3*Q4 +Covar1 +Block, dropone=TRUE, get.ests=TRUE, pheno.col=1, covar=covarAB)       
summary(out.fqB2.a)

result.a <- summary(out.fqB2.a)

rqtl2.a <-refineqtl(RIL2, qtl=qtl2.a, formula=y~Q1+Q2+Q3+Q4+Q5+Q6 +Q1*Q2 + Q3*Q4 +Covar1+Block, verbose=FALSE, pheno.col=1, covar=covarAB)
rqtl2.a
plot(rqtl2.a)
summary(rqtl2.a)

out.fqB3.a <- fitqtl(RIL2, qtl=rqtl2.a, formula=y~Q1+Q2+Q3+Q4+Q5+Q6 +Q1*Q2 + Q3*Q4 +Covar1+Block, dropone=TRUE, get.ests=TRUE, pheno.col=1, covar=covarAB)      
summary(out.fqB3.a)

Everything seems to go OK, except when I run the refineqtl command. There I get twice the error.  I do still get the estimates end results from the summary(out.fqB3.a) command, but not sure how reliable they still are after these errors?
Warning messages:
1: In min(diff(a)) : no non-missing arguments to min; returning Inf
2: In min(diff(a)) : no non-missing arguments to min; returning Inf

For the stepwiseQTL I use the following script:
As you can see I do the stepwise both without giving it the qtl model I built above, and with giving it the qtl model.
In both cases I get the same min(diff(a) errors.

stepperm.0 <- PermStep0
print(pen.0 <- calc.penalties(stepperm.0))
calc.penalties(stepperm.0, alpha=c(0.05, 0.20))
outsw.0 <- stepwiseqtl(stunt2, max.qtl=10, penalties=pen.0, verbose=TRUE, pheno.col=1)
outsw.0qtl <- stepwiseqtl(stunt2, qtl=rqtl2.0, max.qtl=10, penalties=pen.0, verbose=TRUE, pheno.col=1)
summary(outsw.0)

The output from the stepwiseQTL then looks like this:
 -Initial scan
initial lod:  2.827926 
    no.qtl =  1   pLOD = -0.02451691   formula: y ~ Q1 
 -Step 1 
 ---Scanning for additive qtl
        plod = -0.7469474 
 ---Scanning for QTL interacting with Q1
        plod = -0.7027394 
 ---Refining positions
    no.qtl =  2   pLOD = -0.7027394   formula: y ~ Q1 + Q2 + Q1:Q2 
 -Step 2 
 ---Scanning for additive qtl
        plod = -0.8842487 
 ---Scanning for QTL interacting with Q1
        plod = -3.9645 
 ---Scanning for QTL interacting with Q2
        plod = -3.536162 
 ---Look for additional interactions
 ---Refining positions
    no.qtl =  3   pLOD = -0.8842487   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 
 -Step 3 
 ---Scanning for additive qtl
        plod = -0.9683736 
 ---Scanning for QTL interacting with Q1
        plod = -4.418579 
 ---Scanning for QTL interacting with Q2
        plod = -3.000589 
 ---Scanning for QTL interacting with Q3
        plod = -1.843037 
 ---Look for additional interactions
        plod = -4.548217 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  4   pLOD = -0.8705013   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 
 -Step 4 
 ---Scanning for additive qtl
        plod = -1.536615 
 ---Scanning for QTL interacting with Q1
        plod = -3.438925 
 ---Scanning for QTL interacting with Q2
        plod = -3.742803 
 ---Scanning for QTL interacting with Q3
        plod = -2.029572 
 ---Scanning for QTL interacting with Q4
        plod = -2.143622 
 ---Look for additional interactions
        plod = -3.4315 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  5   pLOD = -1.528168   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 
 -Step 5 
 ---Scanning for additive qtl
        plod = -2.205181 
 ---Scanning for QTL interacting with Q1
        plod = -4.196796 
 ---Scanning for QTL interacting with Q2
        plod = -4.487766 
 ---Scanning for QTL interacting with Q3
        plod = -2.627724 
 ---Scanning for QTL interacting with Q4
        plod = -2.986355 
 ---Scanning for QTL interacting with Q5
        plod = -3.173026 
 ---Look for additional interactions
        plod = -4.233457 
 ---Refining positions
    no.qtl =  6   pLOD = -2.205181   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 
 -Step 6 
 ---Scanning for additive qtl
        plod = -2.210587 
 ---Scanning for QTL interacting with Q1
        plod = -5.279551 
 ---Scanning for QTL interacting with Q2
        plod = -5.481877 
 ---Scanning for QTL interacting with Q3
        plod = -2.606345 
 ---Scanning for QTL interacting with Q4
        plod = -3.372474 
 ---Scanning for QTL interacting with Q5
        plod = -3.319364 
 ---Scanning for QTL interacting with Q6
        plod = -2.43476 
 ---Look for additional interactions
        plod = -4.891239 
 ---Refining positions
    no.qtl =  7   pLOD = -2.210587   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 + Q7 
 -Step 7 
 ---Scanning for additive qtl
        plod = -3.011286 
 ---Scanning for QTL interacting with Q1
        plod = -5.882255 
 ---Scanning for QTL interacting with Q2
        plod = -6.259255 
 ---Scanning for QTL interacting with Q3
        plod = -3.850898 
 ---Scanning for QTL interacting with Q4
        plod = -3.32442 
 ---Scanning for QTL interacting with Q5
        plod = -4.061947 
 ---Scanning for QTL interacting with Q6
        plod = -3.870089 
 ---Scanning for QTL interacting with Q7
        plod = -3.236299 
 ---Look for additional interactions
        plod = -2.43476 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  7   pLOD = -2.393031   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q6:Q7 
 -Step 8 
 ---Scanning for additive qtl
        plod = -3.384163 
 ---Scanning for QTL interacting with Q1
        plod = -6.58248 
 ---Scanning for QTL interacting with Q2
        plod = -6.881917 
 ---Scanning for QTL interacting with Q3
        plod = -3.774467 
 ---Scanning for QTL interacting with Q4
        plod = -3.79439 
 ---Scanning for QTL interacting with Q5
        plod = -4.45029 
 ---Scanning for QTL interacting with Q6
        plod = -6.782717 
 ---Scanning for QTL interacting with Q7
        plod = -6.343433 
 ---Look for additional interactions
        plod = -5.314793 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  8   pLOD = -3.370052   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q6:Q7 + Q8 
 -Step 9 
 ---Scanning for additive qtl
        plod = -4.526974 
 ---Scanning for QTL interacting with Q1
        plod = -8.207325 
 ---Scanning for QTL interacting with Q2
        plod = -7.945324 
 ---Scanning for QTL interacting with Q3
        plod = -5.24144 
 ---Scanning for QTL interacting with Q4
        plod = -5.849364 
 ---Scanning for QTL interacting with Q5
        plod = -5.465123 
 ---Scanning for QTL interacting with Q6
        plod = -7.640773 
 ---Scanning for QTL interacting with Q7
        plod = -7.445907 
 ---Scanning for QTL interacting with Q8
        plod = -4.67991 
 ---Look for additional interactions
        plod = -3.729001 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  8   pLOD = -3.661782   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q6:Q7 + Q8 + Q4:Q8 
 -Step 10 
 ---Scanning for additive qtl
        plod = -4.873168 
 ---Scanning for QTL interacting with Q1
        plod = -8.180601 
 ---Scanning for QTL interacting with Q2
        plod = -8.263317 
 ---Scanning for QTL interacting with Q3
        plod = -5.238266 
 ---Scanning for QTL interacting with Q4
        plod = -8.317964 
 ---Scanning for QTL interacting with Q5
        plod = -5.953323 
 ---Scanning for QTL interacting with Q6
        plod = -8.145952 
 ---Scanning for QTL interacting with Q7
        plod = -7.804813 
 ---Scanning for QTL interacting with Q8
        plod = -7.050537 
 ---Look for additional interactions
        plod = -8.659158 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  9   pLOD = -4.469088   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q6:Q7 + Q8 + Q4:Q8 + Q9 
 -Step 11 
 ---Scanning for additive qtl
        plod = -5.49663 
 ---Scanning for QTL interacting with Q1
        plod = -9.108128 
 ---Scanning for QTL interacting with Q2
        plod = -9.09262 
 ---Scanning for QTL interacting with Q3
        plod = -6.356448 
 ---Scanning for QTL interacting with Q4
        plod = -8.828459 
 ---Scanning for QTL interacting with Q5
        plod = -6.570825 
 ---Scanning for QTL interacting with Q6
        plod = -8.600178 
 ---Scanning for QTL interacting with Q7
        plod = -8.184486 
 ---Scanning for QTL interacting with Q8
        plod = -7.973686 
 ---Scanning for QTL interacting with Q9
        plod = -6.521231 
 ---Look for additional interactions
        plod = -9.21341 
 ---Refining positions
 ---  Moved a bit
    no.qtl =  10   pLOD = -5.327537   formula: y ~ Q1 + Q2 + Q1:Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q6:Q7 + Q8 + Q4:Q8 + Q9 + Q10 
 -Starting backward deletion
 ---Dropping Q9 
    no.qtl =  9   pLOD = -4.20007   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q8 + Q9 + Q1:Q2 + Q6:Q7 + Q4:Q8 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q6:Q7 
    no.qtl =  9   pLOD = -3.798451   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q8 + Q9 + Q1:Q2 + Q4:Q8 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q3 
    no.qtl =  8   pLOD = -3.119903   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q8 + Q1:Q2 + Q3:Q7 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q3:Q7 
    no.qtl =  8   pLOD = -2.705104   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q8 + Q1:Q2 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q5 
    no.qtl =  7   pLOD = -2.223793   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q7 + Q1:Q2 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q5 
    no.qtl =  6   pLOD = -1.703857   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q1:Q2 
 ---Refining positions
 ---Dropping Q3 
    no.qtl =  5   pLOD = -1.106742   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q1:Q2 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q4 
    no.qtl =  4   pLOD = -1.222707   formula: y ~ Q1 + Q2 + Q3 + Q4 + Q1:Q2 
 ---Refining positions
 ---  Moved a bit
 ---Dropping Q3 
    no.qtl =  3   pLOD = -0.8842487   formula: y ~ Q1 + Q2 + Q3 + Q1:Q2 
 ---Refining positions
 ---Dropping Q3 
    no.qtl =  2   pLOD = -0.7027394   formula: y ~ Q1 + Q2 + Q1:Q2 
 ---Refining positions
 ---Dropping Q2 
    no.qtl =  1   pLOD = -0.02451691   formula: y ~ Q1 
 ---Refining positions

After that, I get 44 warnings, all the "same min(diff(a): no non-missing arguments to min; returning Inf" message.

I also did the stepwise analysis again with my covariates included, and still get the same errors.


Does this information help?
Thanks!

Jasper

Op maandag 16 september 2019 22:02:48 UTC+2 schreef Karl Broman:

Karl Broman

unread,
Sep 18, 2019, 11:45:20 AM9/18/19
to R/qtl discussion
Do you maybe have a chromosome with a single marker?

karl

On Tuesday, September 17, 2019 at 11:18:05 AM UTC-5, jspv2 v wrote:
Hi Karl,

As in my PhD, I'm working together with several companies, the data is confidential, so not sure if I can share the real data here.
But I can share with you the code/script I used:
 [clip]

jspv2 v

unread,
Sep 18, 2019, 2:38:06 PM9/18/19
to R/qtl discussion
Hi Karl,

Yes, on 2 chromosomes there was such a large gap between markers that these chromosomes split up into 2 linkage groups. (If I force them together, I get a massive gap between them) So there are 2 cases where 1 markers is seen as a whole chromosome. 

I now removed these 2 markers, and the warnings dissapeared!
So this analysis don't work when chromosomes have only a single marker?

Thanks!

Best Wishes,
Jasper

Op woensdag 18 september 2019 17:45:20 UTC+2 schreef Karl Broman:

Karl Broman

unread,
Sep 18, 2019, 2:43:37 PM9/18/19
to rqtl...@googlegroups.com
I think the analysis is okay, and that the warning is bug. And I think I've got it fixed.

   https://github.com/kbroman/qtl/commit/88a81df6908cf2ccfcec235fbc4a628684732007#diff-3f63f45841e7977e6eb9b96dff88f692L141-R144

karl

--
You received this message because you are subscribed to the Google Groups "R/qtl discussion" group.
To unsubscribe from this group and stop receiving emails from it, send an email to rqtl-disc+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/rqtl-disc/caeebe4a-7697-45c2-8bec-2ef0cc998157%40googlegroups.com.

jspv2 v

unread,
Sep 23, 2019, 4:47:59 AM9/23/19
to R/qtl discussion
Hi Karl,

Indeed, it seems to work now!

Thanks so much!

Best Wishes,
Jasper

Op woensdag 18 september 2019 20:43:37 UTC+2 schreef Karl Broman:
I think the analysis is okay, and that the warning is bug. And I think I've got it fixed.

   https://github.com/kbroman/qtl/commit/88a81df6908cf2ccfcec235fbc4a628684732007#diff-3f63f45841e7977e6eb9b96dff88f692L141-R144

karl


On 9/18/19 1:38 PM, 'jspv2 v' via R/qtl discussion wrote:
Hi Karl,

Yes, on 2 chromosomes there was such a large gap between markers that these chromosomes split up into 2 linkage groups. (If I force them together, I get a massive gap between them) So there are 2 cases where 1 markers is seen as a whole chromosome. 

I now removed these 2 markers, and the warnings dissapeared!
So this analysis don't work when chromosomes have only a single marker?

Thanks!

Best Wishes,
Jasper

Op woensdag 18 september 2019 17:45:20 UTC+2 schreef Karl Broman:
Do you maybe have a chromosome with a single marker?

karl

On Tuesday, September 17, 2019 at 11:18:05 AM UTC-5, jspv2 v wrote:
Hi Karl,

As in my PhD, I'm working together with several companies, the data is confidential, so not sure if I can share the real data here.
But I can share with you the code/script I used:
 [clip]
--
You received this message because you are subscribed to the Google Groups "R/qtl discussion" group.
To unsubscribe from this group and stop receiving emails from it, send an email to rqtl...@googlegroups.com.

jspv2 v

unread,
Oct 17, 2019, 5:26:51 AM10/17/19
to R/qtl discussion
Hi Karl,

Now that I have performed the stepwiseQTL analysis, I get some contradictory results.

So in the scanone, I don't find significant QTLs. In the scantwo, i find some pairs of interacting QTLs, that have LODs that are slightly higher than the significance tresholds.  I combined these QTLs in a Multiple QTL model and then used refineQTL and fitQTL to get the results for the complete model.
Here I get several QTLs that are clearly over the significance threshold, so that are clearly significant.
Strangely when I then do the stepwiseQTL, no significant QTLs model comes out at all.

See below the results from the refineQTL & fitQTL, and the significance threshold (from scantwo permutation test with method="em"). (I had to mask the positions of the QTLs with letters due to the confidentiality of my PhD project)

Which result should I trust now, the results from the stepwiseQTL, or from the refineQTL & fitQTL? Also, what I dont completely understand yet, how do the LOD-scores from this multiple QTL model compare to the LOD-scores from the simple scanone analysis? As with scanone I dont get any really significant peaks? How can there be such differences between the different methods?

Thanks!

Best Wishes,
Jasper 

Perm (1008 permutations)
    full  fv1  int  add  av1  one
20% 5.33 4.21 3.28 4.15 2.85 2.16
5%  6.10 4.78 3.90 4.89 3.38 2.76
1%  6.78 5.33 4.46 5.52 4.00 3.47

Full model result
----------------------------------  
Model formula: y ~ Q1 + Q2 + Q3 + Q4 + Q5 + Q6 + Q1:Q2 + Q3:Q4 

       df       SS        MS      LOD     %var Pvalue(Chi2)    Pvalue(F)
Model   8 10948.65 1368.5808 16.56774 31.45702 2.711165e-13 7.853718e-13
Error 193 23856.45  123.6086                                            
Total 201 34805.10                                                      


Drop one QTL at a time ANOVA table: 
----------------------------------  
              df Type III SS   LOD   %var F value Pvalue(Chi2) Pvalue(F)    
x@y            2        4277 7.233 12.287  17.299        0.000  1.23e-07 ***
x@y            2        3317 5.711  9.531  13.418        0.000  3.50e-06 ***
x@y            2        2604 4.544  7.481  10.532        0.000  4.56e-05 ***
x@y            2        2894 5.022  8.315  11.706        0.000  1.59e-05 ***
x@y            1        1468 2.619  4.217  11.875        0.001  0.000698 ***
x@y            1        1142 2.050  3.280   9.236        0.002  0.002702 ** 
x@y:a@b        1        3197 5.516  9.184  25.860        0.000  8.65e-07 ***
x@y:a@b        1        2561 4.473  7.359  20.722        0.000  9.39e-06 ***

Op maandag 23 september 2019 10:47:59 UTC+2 schreef jspv2 v:

Karl Broman

unread,
Oct 17, 2019, 7:18:05 AM10/17/19
to rqtl...@googlegroups.com
Here you are including interactions, whereas in the scanone() analysis you were not. I cannot comment on which is more trustworthy. You should definitely look at plots of phenotype x genotype. 

karl

On Oct 17, 2019, at 4:26 AM, 'jspv2 v' via R/qtl discussion <rqtl...@googlegroups.com> wrote:


To unsubscribe from this group and stop receiving emails from it, send an email to rqtl-disc+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/rqtl-disc/7046f495-ae13-4e3b-9ee0-ef391fb516ce%40googlegroups.com.

jspv2 v

unread,
Oct 17, 2019, 7:58:50 AM10/17/19
to R/qtl discussion
Okay, thanks!

What kind of plots do you mean? Like plots with the average phenotype of all allele combinations of the QTLs? I' could try that together with some statistics, to see if I see significant differences between the different allele combinations?

Or do you mean another kind of plot?

Best Wishes,
Jasper


Op donderdag 17 oktober 2019 13:18:05 UTC+2 schreef Karl Broman:

Karl Broman

unread,
Oct 17, 2019, 8:34:08 AM10/17/19
to rqtl...@googlegroups.com
I was thinking plotPXG() with individual data points, considering two QTL at a time.

karl

On Oct 17, 2019, at 6:58 AM, 'jspv2 v' via R/qtl discussion <rqtl...@googlegroups.com> wrote:


To unsubscribe from this group and stop receiving emails from it, send an email to rqtl-disc+...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/rqtl-disc/f1a08114-722c-44d0-bca5-05db2cc52598%40googlegroups.com.
Reply all
Reply to author
Forward
0 new messages