ASMAP error in lods subscript out of bounds
27 views
Skip to first unread message
Siri B
May 30, 2024, 6:47:53 AMMay 30
to R/qtl discussion
Hi,
#remove markers with the most extreme segregation distortion patterns with the bonferroni criterion
> mm <- statMark(AN_cross_new, stat.type = "marker")$marker$AB
> AN_cross_sd_new <- drop.markers(AN_cross_new, c(markernames(AN_cross_new)[mm > 0.98], markernames(AN_cross_new)[mm < 0.2]))
#remove markers with less than 90% present
> AN_cross_pullseg1 <- pullCross(AN_cross_sd_new, type = "missing", pars = list(miss.thresh = 0.1))
#remove colocated markers (NEW)
> AN_cross_pullseg2 <- pullCross(AN_cross_pullseg1, type = "co.located") #Not in Marie's script
#pull out distorted markers with seg.thresh = 1e-10
> AN_cross_pullseg_3<-pullCross(AN_cross_pullseg2, type = "seg.distortion", pars = list(seg.thresh = 1e-10))
#AN_cross_pullseg$seg.distortion$table
> summary(AN_cross_pullseg_3)
BC(0)F(2) cross
No. individuals: 332
No. phenotypes: 6
Percent phenotyped: 100 99.7 99.7 99.7 99.7 99.7
No. chromosomes: 8
Autosomes: 1 112 223 334 445 556 667 762
Total markers: 1351
No. markers: 234 128 192 170 141 194 161 131
Percent genotyped: 94.1
Genotypes (%): AA:26.7 AB:48.9 BB:24.4 not BB:0.0 not AA:0.0
Warning message:
In summary.cross(AN_cross_pullseg_3) :
Some chromosomes > 1000 cM in length; there may be a problem with the genetic map.
(Perhaps it is in basepairs?)
#Making new map without seg. dist markers, co-located markers, and markers with a lot of missing data
> AN_cross_filtered_new <- mstmap.cross(AN_cross_pullseg_3, bychr = TRUE, dist.fun = "kosambi", objective.fun="COUNT", anchor = TRUE, id = "ID", trace = FALSE, mvest.bc=TRUE, return.imputed=TRUE, detectBadData=TRUE,p.value=3)
> summary(AN_cross_filtered_new)
BC(0)F(2) cross
No. individuals: 332
No. phenotypes: 6
Percent phenotyped: 100 99.7 99.7 99.7 99.7 99.7
No. chromosomes: 8
Autosomes: 1 112 223 334 445 556 667 762
Total markers: 1351
No. markers: 234 128 192 170 141 194 161 131
Percent genotyped: 94.1
Genotypes (%): AA:26.7 AB:48.9 BB:24.4 not BB:0.0 not AA:0.0
Warning message:
In summary.cross(AN_cross_filtered_new) :
Some markers at the same position on chr 1,112,223,334,445,556,667,762; use jittermap().
#summary
> summary.map(AN_cross_filtered_new)
n.mar length ave.spacing max.spacing
1 234 863.1 3.7 47.5
112 128 524.7 4.1 83.0
223 192 799.5 4.2 80.3
334 170 658.7 3.9 78.2
445 141 554.2 4.0 91.9
556 194 644.8 3.3 91.5
667 161 611.9 3.8 26.7
762 131 663.2 5.1 67.4
overall 1351 5320.1 4.0 91.9
We have been stuck on a problem with pushcross() in asmap for a while. (I know asmap is tangential to r/qtl, but maybe someone on this forum will be able to help).
Info: We have a lot of seggregation distortion among our markers, which is expected because the study is on hybrid incompatibilities (probably selection against heterozygotes). We therefore want to remove distorted markers (and markers with a lot of missing data + co-located markers) from the map before rebuilding it and pushing the removed markers back in. However, we always end up with the same error message when pushing the markers back into the new map (Error in lods[newmar, newmar] : subscript out of bounds).
There are also a couple of warning messages regarding length + markers at the same position upstream. We are uncertain if we should ignore these or if they might be relevant to our problem.
Any idea what we are doing wrong here?
See code below (I removed some of the data checking commands like plotting for simplicity).
Best wishes,
Siri
Excerpt of our code:
> AN_cross_new<-convert2bcsft(cross1, BC.gen = 0, F.gen = 2)
> summary(AN_cross_new)
BC(0)F(2) cross
No. individuals: 332
No. phenotypes: 6
Percent phenotyped: 100 99.7 99.7 99.7 99.7 99.7
No. chromosomes: 8
Autosomes: 1 112 223 334 445 556 667 762
Total markers: 5437
No. markers: 819 548 675 735 519 842 652 647
Percent genotyped: 91.5
Genotypes (%): AA:36.5 AB:42.4 BB:21.1 not BB:0.0 not AA:0.0
Warning message:
In summary.cross(AN_cross_new) :
Some chromosomes > 1000 cM in length; there may be a problem with the genetic map.
(Perhaps it is in basepairs?)
> summary(AN_cross_new)
BC(0)F(2) cross
No. individuals: 332
No. phenotypes: 6
Percent phenotyped: 100 99.7 99.7 99.7 99.7 99.7
No. chromosomes: 8
Autosomes: 1 112 223 334 445 556 667 762
Total markers: 5437
No. markers: 819 548 675 735 519 842 652 647
Percent genotyped: 91.5
Genotypes (%): AA:36.5 AB:42.4 BB:21.1 not BB:0.0 not AA:0.0
Warning message:
In summary.cross(AN_cross_new) :
Some chromosomes > 1000 cM in length; there may be a problem with the genetic map.
(Perhaps it is in basepairs?)
#remove markers with the most extreme segregation distortion patterns with the bonferroni criterion
> mm <- statMark(AN_cross_new, stat.type = "marker")$marker$AB
> AN_cross_sd_new <- drop.markers(AN_cross_new, c(markernames(AN_cross_new)[mm > 0.98], markernames(AN_cross_new)[mm < 0.2]))
#remove markers with less than 90% present
> AN_cross_pullseg1 <- pullCross(AN_cross_sd_new, type = "missing", pars = list(miss.thresh = 0.1))
#remove colocated markers (NEW)
> AN_cross_pullseg2 <- pullCross(AN_cross_pullseg1, type = "co.located") #Not in Marie's script
#pull out distorted markers with seg.thresh = 1e-10
> AN_cross_pullseg_3<-pullCross(AN_cross_pullseg2, type = "seg.distortion", pars = list(seg.thresh = 1e-10))
#AN_cross_pullseg$seg.distortion$table
> summary(AN_cross_pullseg_3)
BC(0)F(2) cross
No. individuals: 332
No. phenotypes: 6
Percent phenotyped: 100 99.7 99.7 99.7 99.7 99.7
No. chromosomes: 8
Autosomes: 1 112 223 334 445 556 667 762
Total markers: 1351
No. markers: 234 128 192 170 141 194 161 131
Percent genotyped: 94.1
Genotypes (%): AA:26.7 AB:48.9 BB:24.4 not BB:0.0 not AA:0.0
Warning message:
In summary.cross(AN_cross_pullseg_3) :
Some chromosomes > 1000 cM in length; there may be a problem with the genetic map.
(Perhaps it is in basepairs?)
#Making new map without seg. dist markers, co-located markers, and markers with a lot of missing data
> AN_cross_filtered_new <- mstmap.cross(AN_cross_pullseg_3, bychr = TRUE, dist.fun = "kosambi", objective.fun="COUNT", anchor = TRUE, id = "ID", trace = FALSE, mvest.bc=TRUE, return.imputed=TRUE, detectBadData=TRUE,p.value=3)
> summary(AN_cross_filtered_new)
BC(0)F(2) cross
No. individuals: 332
No. phenotypes: 6
Percent phenotyped: 100 99.7 99.7 99.7 99.7 99.7
No. chromosomes: 8
Autosomes: 1 112 223 334 445 556 667 762
Total markers: 1351
No. markers: 234 128 192 170 141 194 161 131
Percent genotyped: 94.1
Genotypes (%): AA:26.7 AB:48.9 BB:24.4 not BB:0.0 not AA:0.0
Warning message:
In summary.cross(AN_cross_filtered_new) :
Some markers at the same position on chr 1,112,223,334,445,556,667,762; use jittermap().
#summary
> summary.map(AN_cross_filtered_new)
n.mar length ave.spacing max.spacing
1 234 863.1 3.7 47.5
112 128 524.7 4.1 83.0
223 192 799.5 4.2 80.3
334 170 658.7 3.9 78.2
445 141 554.2 4.0 91.9
556 194 644.8 3.3 91.5
667 161 611.9 3.8 26.7
762 131 663.2 5.1 67.4
overall 1351 5320.1 4.0 91.9
# Check dimensions and markers
> print(dim(AN_cross_filtered_new))
> print(dim(AN_cross_filtered_new))
NULL
> print(totmar(AN_cross_filtered_new))
> print(totmar(AN_cross_filtered_new))
[1] 1351
# Recalculate LOD scores and recombination fractions
> AN_cross_filtered_new <- est.rf(AN_cross_filtered_new)
#push back distorted markers with seg.thresh = 1e-20
#create map for QTL analysis
> AN_cross_f2_no_pushseg_new<-pushCross(AN_cross_filtered_new, type = "seg.distortion", pars = list(seg.thresh = 1e-20))
Error in lods[newmar, newmar] : subscript out of bounds
# Recalculate LOD scores and recombination fractions
> AN_cross_filtered_new <- est.rf(AN_cross_filtered_new)
#push back distorted markers with seg.thresh = 1e-20
#create map for QTL analysis
> AN_cross_f2_no_pushseg_new<-pushCross(AN_cross_filtered_new, type = "seg.distortion", pars = list(seg.thresh = 1e-20))
Error in lods[newmar, newmar] : subscript out of bounds
Reply all
Reply to author
Forward
0 new messages