Skip to content

Commit 275d978

Browse files
gregorgorjancgregorgorjanc
committed
Merge remote-tracking branch 'upstream/devel' into ts
2 parents 22c13a9 + 5fc269f commit 275d978

File tree

5 files changed

+67
-28
lines changed

5 files changed

+67
-28
lines changed

‎NEWS.md‎

Lines changed: 5 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,5 @@
1+
*fixed bug in `getGvE` with multiple traits.
2+
13
# AlphaSimR 2.1.0
24

35
*changed R6 and methods from Depends to Imports to match current best practices for R packages
@@ -209,7 +211,7 @@
209211
*removed lazyData field in DESCRIPTION
210212

211213
# AlphaSimR 1.0.0
212-
214+
213215
*AlphaSimR manuscript has been published in G3 (citation added)
214216

215217
*changed to a Gamma Sprinkling model for crossovers, default is still a Gamma model
@@ -280,7 +282,7 @@
280282

281283
*the `c` function now merges individuals for MapPop objects (was chromosomes before)
282284

283-
*the `cChr` function new merges chromosomes for MapPop objects
285+
*the `cChr` function new merges chromosomes for MapPop objects
284286

285287
*fixed broken SimParam_addStructuredSnpChip
286288

@@ -447,4 +449,4 @@
447449
*Added `runMacs2` as a wrapper for `runMacs`
448450

449451
*Fixed error when using H2 in SimParam_setVarE
450-
452+

‎R/GS.R‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -14,8 +14,8 @@ convertTraitsToNames = function(traits, simParam){
1414
if(is.character(traits)){
1515
# Suspect trait is a name
1616
take = match(traits, simParam$traitNames)
17-
if(is.na(take)){
18-
stop("'",traits,"' did not match any trait names")
17+
if(any(is.na(take))){
18+
stop("'",traits[is.na(take)],"' did not match any trait names")
1919
}
2020
traits = take
2121
}else if(is.function(traits)){

‎R/selection.R‎

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -37,8 +37,8 @@ getResponse = function(pop,trait,use,simParam=NULL,...){
3737
}else{ # trait is not a function, so must be numeric or character
3838
if(is.character(trait)){ # Suspect trait is a name
3939
take = match(trait, simParam$traitNames)
40-
if(is.na(take)){
41-
stop("'",trait,"' did not match any trait names")
40+
if(any(is.na(take))){
41+
stop("'",trait[is.na(take)],"' did not match any trait names")
4242
}
4343
trait = take
4444
}

‎src/getGv.cpp‎

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -351,11 +351,11 @@ arma::field<arma::vec> getGvIndexE(const Rcpp::S4& trait,
351351
d(E(i,0))*xd(genoMat(j,qtlIndex(E(i,0)))) +
352352
a(E(i,1))*xa(genoMat(j,qtlIndex(E(i,1)))) +
353353
d(E(i,1))*xd(genoMat(j,qtlIndex(E(i,1)))) +
354-
E(i,2)*xa(genoMat(j,(E(i,0))))*xa(genoMat(j,qtlIndex(E(i,1))));
354+
E(i,2)*xa(genoMat(j,qtlIndex(E(i,0))))*xa(genoMat(j,qtlIndex(E(i,1))));
355355
}else{
356356
gv(j,tid) += a(E(i,0))*xa(genoMat(j,qtlIndex(E(i,0)))) +
357-
a(E(i,1))*xa(genoMat(j,(E(i,1)))) +
358-
E(i,2)*xa(genoMat(j,(E(i,0))))*xa(genoMat(j,qtlIndex(E(i,1))));
357+
a(E(i,1))*xa(genoMat(j,qtlIndex(E(i,1)))) +
358+
E(i,2)*xa(genoMat(j,qtlIndex(E(i,0))))*xa(genoMat(j,qtlIndex(E(i,1))));
359359
}
360360
if(hasGxe){
361361
gxe(j,tid) += g(E(i,0))*xa(genoMat(j,qtlIndex(E(i,0)))) +

‎tests/testthat/test-phenotypes.R‎

Lines changed: 55 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -1,34 +1,44 @@
11
context("phenotypes")
22

3-
test_that("asCategorical_converts_correctly",{
3+
test_that("asCategorical_converts_correctly", {
44
cont = matrix(data = 0, nrow = 7, ncol = 3)
55
cont[, 1] = c(-3, -2, -1, 0, 1, 2, 3)
66
cont[, 2] = c(-3, -2, -1, 0, 1, 2, 3)
77
cont[, 3] = c(-3, -2, -1, 0, 1, 2, 3)
88

9-
expect_equal(asCategorical(x = cont[, 1]),
10-
matrix(c(1, 1, 1, 2, 2, 2, 2)))
11-
expect_equal(asCategorical(x = cont[, 1], threshold = c(-1, 0, 1)),
12-
matrix(c(NA, NA, 1, 2, 2, NA, NA)))
13-
expect_equal(asCategorical(x = cont[, 1], threshold = c(-Inf, -1, 0, 1, Inf)),
14-
matrix(c(1, 1, 2, 3, 4, 4, 4)))
9+
expect_equal(asCategorical(x = cont[, 1]), matrix(c(1, 1, 1, 2, 2, 2, 2)))
10+
expect_equal(
11+
asCategorical(x = cont[, 1], threshold = c(-1, 0, 1)),
12+
matrix(c(NA, NA, 1, 2, 2, NA, NA))
13+
)
14+
expect_equal(
15+
asCategorical(x = cont[, 1], threshold = c(-Inf, -1, 0, 1, Inf)),
16+
matrix(c(1, 1, 2, 3, 4, 4, 4))
17+
)
1518

1619
expect_warning(asCategorical(x = cont[, 1], p = 0.5))
17-
expect_equal(suppressWarnings(asCategorical(x = cont[, 1], p = 0.5)),
18-
asCategorical(x = cont[, 1], p = c(0.5, 0.5)))
20+
expect_equal(
21+
suppressWarnings(asCategorical(x = cont[, 1], p = 0.5)),
22+
asCategorical(x = cont[, 1], p = c(0.5, 0.5))
23+
)
1924

2025
trtMean = apply(X = cont, MARGIN = 2, FUN = mean)
2126
trtVar = apply(X = cont, MARGIN = 2, FUN = var)
22-
expect_equal(asCategorical(x = cont[, 1], p = c(0.5, 0.5), var = trtVar[1]),
23-
matrix(c(1, 1, 1, 2, 2, 2, 2)))
24-
expect_equal(asCategorical(x = cont[, 1], p = c(2/7, 1/7, 1/7, 3/7), var = trtVar[1]),
25-
matrix(c(1, 1, 2, 3, 4, 4, 4)))
27+
expect_equal(
28+
asCategorical(x = cont[, 1], p = c(0.5, 0.5), var = trtVar[1]),
29+
matrix(c(1, 1, 1, 2, 2, 2, 2))
30+
)
31+
expect_equal(
32+
asCategorical(
33+
x = cont[, 1],
34+
p = c(2 / 7, 1 / 7, 1 / 7, 3 / 7),
35+
var = trtVar[1]
36+
),
37+
matrix(c(1, 1, 2, 3, 4, 4, 4))
38+
)
2639

2740
expect_error(asCategorical(x = cont))
28-
cont2 = asCategorical(x = cont,
29-
threshold = list(NULL,
30-
c(-Inf, 0, Inf),
31-
NULL))
41+
cont2 = asCategorical(x = cont, threshold = list(NULL, c(-Inf, 0, Inf), NULL))
3242
cont2Exp = cont
3343
cont2Exp[, 2] = c(1, 1, 1, 2, 2, 2, 2)
3444
expect_equal(cont2, cont2Exp)
@@ -37,8 +47,35 @@ test_that("asCategorical_converts_correctly",{
3747
pList = list(NULL, c(0.5, 0.5), NULL)
3848
expect_error(asCategorical(x = cont, p = pList))
3949
expect_error(asCategorical(x = cont, p = pList, mean = trtMean))
40-
cont2 = asCategorical(x = cont, p = pList, mean = trtMean, var=trtVar)
50+
cont2 = asCategorical(x = cont, p = pList, mean = trtMean, var = trtVar)
4151
cont2Exp = cont
4252
cont2Exp[, 2] = c(1, 1, 1, 2, 2, 2, 2)
4353
expect_equal(cont2, cont2Exp)
4454
})
55+
56+
test_that("pop@gv and genParam(pop)@gv match", {
57+
founderPop = quickHaplo(nInd = 10, nChr = 1, segSites = 10, ploidy = 1)
58+
SP = SimParam$new(founderPop)
59+
SP$addTraitA(nQtlPerChr = 10, mean = 0, var = 1, name = "addTraitA_allQTLs")
60+
SP$addTraitAE(
61+
nQtlPerChr = 10,
62+
relAA = 0,
63+
mean = 0,
64+
var = 1,
65+
useVarA = FALSE,
66+
name = "addTraitAE_allQTLs"
67+
)
68+
SP$addTraitA(nQtlPerChr = 2, mean = 0, var = 1, name = "addTraitA_2QTLs")
69+
SP$addTraitAE(
70+
nQtlPerChr = 2,
71+
relAA = 0,
72+
mean = 0,
73+
var = 1,
74+
useVarA = FALSE,
75+
name = "addTraitAE_2QTLs"
76+
)
77+
pop = newPop(founderPop, simParam = SP)
78+
diff = pop@gv - genParam(pop, simParam = SP)$gv
79+
test = abs(diff) < .Machine$double.eps^0.5
80+
expect_true(all(test))
81+
})

0 commit comments

Comments
 (0)