singular matricex using cim()

[andbelo]

Dear Karl,

I have been using R/qtl for mapping and found a error that happened
some times. I use RILs with a dense genetic map. Using different
phenotypes with the same mapping population I get the error in some of
the cases, so I imagine that the population and marker data is fine.
Also, they work well for interval mapping and for other mapping using
either qtl_cartographer or R/qtlbim. Interestingly, I rerun the same
commands for a couple of phenotypes and got two distinct results: in
some I did not see the error and all the mapping went ok, but in most
of the cases the error appeared again. I also tried changing the
window value for cim but id did not work either.

It seems to me that the error occurs in the execution of the cim()
function. The proccess is generating singular matrices. The plot
resulting of the cim() analysis is a flat QTL line with very high LOD
score (e.g. ~48) and one or two spikes towards more normal LOD score
values.

I would like to know if you have seen this error before and/or know
how to solve it. I would gladly help in any checking that you think
might be useful to troubleshoot this.

Thank you



Bellow I added the relevant commands and output:

R> library(qtl)
R> d = read.cross('qtlcart', file = 'qtlcart.cro', mapfile =
'qtlcart.map')
R> d1 = calc.genoprob(d, error.prob = 0.001, map.function = 'haldane')

... other mapping commands using interval mapping

R> cim = cim(d1, pheno.col = i, window = 10, method = 'em', error.prob
= 0.0001, map.function = 'haldane')
There were 50 or more warnings (use warnings() to see the first 50)
R> warnings();
Warning messages:
1: In checkcovar(cross, pheno.col, addcovar, intcovar,
perm.strata, ... :
addcovar appears to be over-specified; consider dropping columns.

2: X'X matrix is singular.
3: X'X matrix is singular.

... and the same message goes up to:

49: X'X matrix is singular.
50: X'X matrix is singular.

[Karl Broman]

There are probably multiple markers with identical genotypes, and
the forward selection algorithm in cim() [which is rather roughly
implemented] is getting singularities.

You might try first running calc.genoprob with error.prob=0.01 or
=0.001.
That would allow a bit of fuzziness and might prevent this problem.

karl

[Miguel Macias]

Hello Dr. Broman,

 

I am having a similar issue when I run the cim() function using the "ehk" method. The code is given below:

 

#### loading the data and calculating genoprob. I am running the analysis on a RIL population by selfing.

Rqtl.data = read.cross("csv", "F:/LAB_WORK/WUENUE/GRxI_QTL_R/", "my_pheno_geno_data.csv", alleles = c("A", "B"), genotypes = c("AA", "BB"))

Rqtl.data <- convert2riself(Rqtl.data)

grxi.genoprob = calc.genoprob(Rqtl.data, step = 1, error.prob = 0.001)

 

#### running cim() on all my traits

phen.col = 2:35

n.phe = length(phen.col)

 ehk.cim <- ehk.cim.p1000 <- vector("list", n.phe)

 thr.ehk <- rep(NA, n.phe)

 

for(i in seq(1, length(phen.col), 1)) {

 

    ehk.cim.p1000[[i]] <- cim(grxi.genoprob, pheno.col = phen.col[i], n.marcovar = 10, window = 15, method = "ehk", n.perm = 20)

  thr.ehk[i] <- summary(ehk.cim.p1000[[i]], alpha = 0.05)

  ehk.cim[[i]] <- cim(grxi.genoprob, pheno.col = phen.col[i], n.marcovar = 10, window = 15, method = "ehk")

 

}

 

I get 50 warnings that read:

X'X matrix is singular

 

And in some warning messages I can get:

Didn't converge at position 141

The position that didn't converge changes each time I run the cim() function

 

I also tried running calc.genoprob with error.prob = 0.01 but I still get the same warning messages.

 

Here are some diagnostics I have done:

1)     I do not have any markers with the same allele for all individuals (total 89 individuals)

2)     I do not have any markers in the same genetic position

3)     I do not have markers that are linked in different chromosomes

4)     I do have a pair of individuals that are 90% similar and another pair that are 80% similar genotypically

5)     I have 8 RILs that have between 30 and 40 crossovers

 

Can you please provide further advice on how to troubleshoot this problem?

 

Thank you,

 

Miguel Macias Gonzalez

[Karl Broman]

method="ehk" can be a bit unstable, unfortunately. I'd recommend choosing a different method,
such as "hk".

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 post to this group, send email to rqtl...@googlegroups.com.
> Visit this group at http://groups.google.com/group/rqtl-disc.
> For more options, visit https://groups.google.com/d/optout.