non-numeric argument error in boot.ppfis

[jaflury]

Dear Jerome,

I had the same issue with the population name in basic.stats and it was resolved with the updated function, however I get the same error in boot.ppfis (error message "data[dim(data)[1], 1] + 1 : non-numeric argument to binary operator"). Unfortunately I could not figure out what needs to be changed, could you help me there?

Best, Jana

Originally posted by @jaflury in #65 (comment)

[jgx65]

Dear Jana,

Could you post a reproducible example please?

A quick fix is to duplicate your unique population:

#load gtrunchier dataset
data(gtrunchier)
#extract individuals from 1st pop
dum<-gtrunchier[gtrunchier[,1]==1,-2]
#produces an error
boot.ppfis(dum)
#duplicates dum
dum1<-dum
dum1[,1]<-2
dumt<-rbind(dum,dum1)
#works, repeating twice the results
boot.ppfis(dumt)

[jgx65]

@jaflury , have you managed to solve the problem?

[jaflury]

Dear Jerome,
sorry for the late reply. Here's a subset of the data I used. I tried your solution:

subset <- read.table("dat.txt")

dum<-subset[subset[,1]==1,-2]
#produces an error
boot.ppfis(dum)
#duplicates dum
dum1<-dum
dum1[,1]<-2
dumt<-rbind(dum,dum1)
#works, repeating twice the results
boot.ppfis(dumt)

but it didn't work. I know it's coming from the fact that there is only one population, but I don't find a way to fix it. It was fixed in this basic.stat function:

basic.stats<-function (data, diploid = TRUE, digits = 4)
{
if (adegenet::is.genind(data))
data <- genind2hierfstat(data)
if (length(table(data[, 1])) < 2){
data[dim(data)[1] + 1, 1] <- "DumPop"
dum.pop<-TRUE
}
if (dim(data)[2] == 2)
data <- data.frame(data, dummy.loc = data[, 2])
loc.names <- names(data)[-1]
p <- pop.freq(data, diploid)
n <- t(ind.count(data))
if (diploid) {
dum <- getal.b(data[, -1])
Ho <- dum[, , 1] == dum[, , 2]
sHo <- (1 - t(apply(Ho, 2, fun <- function(x) tapply(x,
data[, 1], mean, na.rm = TRUE))))
mHo <- apply(sHo, 1, mean, na.rm = TRUE)
}
else {
sHo <- NA
mHo <- NA
}
sp2 <- lapply(p, fun <- function(x) apply(x, 2, fun2 <- function(x) sum(x^2)))
sp2 <- matrix(unlist(sp2), nrow = dim(data[, -1])[2], byrow = TRUE)
if (diploid) {
Hs <- (1 - sp2 - sHo/2/n)
Hs <- n/(n - 1) * Hs
Fis = 1 - sHo/Hs
}
else {
Hs <- n/(n - 1) * (1 - sp2)
Fis <- NA
}
np <- apply(n, 1, fun <- function(x) sum(!is.na(x)))
mn <- apply(n, 1, fun <- function(x) {
np <- sum(!is.na(x))
np/sum(1/x[!is.na(x)])
})
msp2 <- apply(sp2, 1, mean, na.rm = TRUE)
mp <- lapply(p, fun <- function(x) apply(x, 1, mean, na.rm = TRUE))
mp2 <- unlist(lapply(mp, fun1 <- function(x) sum(x^2)))
if (diploid) {
mHs <- mn/(mn - 1) * (1 - msp2 - mHo/2/mn)
Ht <- 1 - mp2 + mHs/mn/np - mHo/2/mn/np
mFis = 1 - mHo/mHs
}
else {
mHs <- mn/(mn - 1) * (1 - msp2)
Ht <- 1 - mp2 + mHs/mn/np
mFis <- NA
}
Dst <- Ht - mHs
Dstp <- np/(np - 1) * Dst
Htp = mHs + Dstp
Fst = Dst/Ht
Fstp = Dstp/Htp
Dest <- Dstp/(1 - mHs)
res <- data.frame(cbind(mHo, mHs, Ht, Dst, Htp, Dstp, Fst,
Fstp, mFis, Dest))
names(res) <- c("Ho", "Hs", "Ht", "Dst", "Htp", "Dstp",
"Fst", "Fstp", "Fis", "Dest")
if (diploid) {
rownames(sHo) <- loc.names
rownames(Fis) <- loc.names
}
is.na(res) <- do.call(cbind, lapply(res, is.infinite))
overall <- apply(res, 2, mean, na.rm = TRUE)
overall[7] <- overall[4]/overall[3]
overall[8] <- overall[6]/overall[5]
overall[9] <- 1 - overall[1]/overall[2]
overall[10] <- overall[6]/(1 - overall[2])
names(overall) <- names(res)
if (!diploid) {
overall[c(1, 9)] <- NA
}
if(dum.pop){
ToBeRemoved<-which(dimnames(Hs)[[2]]=="DumPop")
n<-n[,-ToBeRemoved,drop=FALSE]
Hs<-Hs[,-ToBeRemoved,drop=FALSE]
if (diploid) sHo<-sHo[,-ToBeRemoved,drop=FALSE]
Fis<-Fis[,-ToBeRemoved,drop=FALSE]
p<-lapply(p,function(x) x[,-ToBeRemoved,drop=FALSE])
}
all.res <- list(n.ind.samp = n, pop.freq = lapply(p, round,
digits), Ho = round(sHo, digits), Hs = round(Hs, digits),
Fis = round(Fis, digits), perloc = round(res, digits),
overall = round(overall, digits))
class(all.res) <- "basic.stats"
all.res
}

and this function works with the dat.txt.
basic.stats(subset)

Thank you very much for your help!

Best, Jana

[jaflury]

and here's the file:
dat.txt

[jgx65]

Dear @jaflury , sorry for my delayed answer. the following works for me:

library(hierfstat)
dat<-read.table("dat.txt",header=T)
tmp<-dat
tmp[,1]<-rep("pop2",dim(tmp)[1])
dat2<-rbind(dat,tmp)
boot.ppfis(dat2)

and I obtain:

$fis.ci
      ll     hl
1 0.0863 0.3339
2 0.0863 0.3339

Let me know if this solves the issue. I will look into a more permanent solution when I find time. Many thanks!