Skip to content

Commit 8ae49d1

Browse files
philchalmersphilchalmers
committed
use simpler logic throughout
1 parent a436ba2 commit 8ae49d1

2 files changed

Lines changed: 59 additions & 131 deletions

File tree

‎R/SingleGroup-methods.R‎

Lines changed: 55 additions & 129 deletions
Original file line numberDiff line numberDiff line change
@@ -1696,6 +1696,7 @@ setMethod(
16961696
)
16971697

16981698
mirt2traditional <- function(x, vcov, nfact){
1699+
xorg <- x
16991700
cls <- class(x)
17001701
opar <- par <- x@par
17011702
which.a <- which(x@par[1L:nfact] != 0)
@@ -1705,54 +1706,23 @@ mirt2traditional <- function(x, vcov, nfact){
17051706
ncat <- x@ncat
17061707
if(cls == 'dich'){
17071708
fns <- vector('list', nfact + 3L)
1708-
fns[[nfact+1L]] <- function(par, index, opar){
1709-
if(index == (nfact + 1L)){
1710-
opar[c(which.a, nfact + 1L)] <- par
1711-
ret <- -opar[nfact + 1L]/opar[which.a]
1712-
}
1713-
ret
1714-
}
1715-
fns[[nfact+2L]] <- function(par, index, opar){
1716-
if(index == nfact + 2L)
1717-
ret <- plogis(par)
1718-
ret
1719-
}
1720-
fns[[nfact+3L]] <- function(par, index, opar){
1721-
if(index == nfact + 3L)
1722-
ret <- plogis(par)
1723-
ret
1724-
}
1709+
fns[[nfact+1L]] <- function(par, index, opar) -opar[2]/par[1]
1710+
fns[[nfact+2L]] <- function(par, index, opar) plogis(par)
1711+
fns[[nfact+3L]] <- function(par, index, opar) plogis(par)
17251712
delta_index <- c(as.list(rep(NA, nfact)),
17261713
list(c(which.a, nfact + 1L),
17271714
nfact + 2L, nfact+3L))
17281715
par[nfact + 1L] <- -par[nfact + 1L]/par[which.a]
17291716
par[nfact + 2L] <- plogis(par[nfact + 2L])
17301717
par[nfact + 3L] <- plogis(par[nfact + 3L])
17311718
names(par) <- c(a.nms, 'b', 'g', 'u')
1719+
index <- delta_index
17321720
} else if(cls == 'fivePL'){
17331721
fns <- vector('list', nfact + 4L)
1734-
fns[[nfact+1L]] <- function(par, index, opar){
1735-
if(index == (nfact + 1L)){
1736-
opar[c(which.a, nfact + 1L)] <- par
1737-
ret <- -opar[nfact + 1L]/opar[which.a]
1738-
}
1739-
ret
1740-
}
1741-
fns[[nfact+2L]] <- function(par, index, opar){
1742-
if(index == nfact + 2L)
1743-
ret <- plogis(par)
1744-
ret
1745-
}
1746-
fns[[nfact+3L]] <- function(par, index, opar){
1747-
if(index == nfact + 3L)
1748-
ret <- plogis(par)
1749-
ret
1750-
}
1751-
fns[[nfact+4L]] <- function(par, index, opar){
1752-
if(index == nfact + 4L)
1753-
ret <- exp(par)
1754-
ret
1755-
}
1722+
fns[[nfact+1L]] <- function(par, index, opar) -opar[2]/par[1]
1723+
fns[[nfact+2L]] <- function(par, index, opar) plogis(par)
1724+
fns[[nfact+3L]] <- function(par, index, opar) plogis(par)
1725+
fns[[nfact+4L]] <- function(par, index, opar) exp(par)
17561726
delta_index <- c(as.list(rep(NA, nfact)),
17571727
list(c(which.a, nfact + 1L),
17581728
nfact + 2L, nfact+3L, nfact+4L))
@@ -1761,17 +1731,11 @@ mirt2traditional <- function(x, vcov, nfact){
17611731
par[nfact + 3L] <- plogis(par[nfact + 3L])
17621732
par[nfact + 4L] <- exp(par[nfact + 4L])
17631733
names(par) <- c(a.nms, 'b', 'g', 'u', 'S')
1734+
index <- delta_index
17641735
} else if(cls == 'graded'){
17651736
fns <- vector('list', ncat + nfact-1L)
1766-
for(i in 2L:ncat - 1L){
1767-
fns[[i + nfact]] <- function(par, index, opar){
1768-
if(index > nfact){
1769-
opar[c(which.a, index)] <- par
1770-
ret <- -opar[index]/opar[which.a]
1771-
}
1772-
ret
1773-
}
1774-
}
1737+
for(i in 2L:ncat - 1L)
1738+
fns[[i + nfact]] <- function(par, index, opar) -par[2]/par[1]
17751739
delta_index <- vector('list', ncat + nfact - 1L)
17761740
for(i in 1:nfact)
17771741
delta_index[[i]] <- NA
@@ -1780,24 +1744,20 @@ mirt2traditional <- function(x, vcov, nfact){
17801744
delta_index[[i+nfact]] <- c(which.a, i+nfact)
17811745
}
17821746
names(par) <- c(a.nms, paste0('b', 2:ncat-1L))
1747+
index <- delta_index
17831748
} else if(cls == 'gpcm'){
17841749
fns <- vector('list', ncat+nfact)
17851750
for(i in 2L:ncat-1L){
17861751
fns[[i+nfact]] <- function(par, index, opar){
1787-
if(index > nfact){
1788-
if(index == (nfact+1))
1789-
opar[c(which.a, ncat + nfact + 2L)] <- par
1790-
else opar[c(which.a, ncat + index + nfact - 1L,
1791-
ncat + index + nfact)] <- par
1792-
par <- opar
1793-
ds <- par[(nfact+1):length(par)]/par[which.a]
1794-
ds <- ds[-seq_len(ncat)]
1795-
newd <- numeric(length(ds)-1L)
1796-
for(i in 2:length(ds))
1797-
newd[i-1L] <- -(ds[i] - ds[i-1L])
1798-
ret <- c(par[1:nfact], newd)
1799-
ret <- ret[index]
1800-
}
1752+
opar[index] <- par
1753+
par <- opar
1754+
ds <- par[(nfact+1):length(par)]/par[which.a]
1755+
ds <- ds[-seq_len(ncat)]
1756+
newd <- numeric(length(ds)-1L)
1757+
for(i in 2:length(ds))
1758+
newd[i-1L] <- -(ds[i] - ds[i-1L])
1759+
ret <- c(par[1:nfact], 1:ncat-1, newd)
1760+
ret <- ret[index[2]]
18011761
ret
18021762
}
18031763
}
@@ -1818,21 +1778,18 @@ mirt2traditional <- function(x, vcov, nfact){
18181778
par <- c(x@par[1:nfact], newd)
18191779
names(par) <- c(a.nms, paste0('b', 1:length(newd)))
18201780
x@est <- x@est[c(1:nfact, (ncat+nfact+2L):length(x@est))]
1781+
index <- delta_index
18211782
} else if(cls == 'nominal'){
1822-
if(nfact > 1L) return(x)
18231783
fns <- vector('list', ncat*2)
1824-
for(i in 2L:length(par)-1L){
1784+
for(i in 1:ncat){
18251785
fns[[i]] <- function(par, index, opar){
1826-
if(index <= floor(length(opar)/2)){
1827-
as <- par[2:(ncat+1)] * par[1]
1828-
as <- as - mean(as)
1829-
ret <- as[index]
1830-
} else {
1831-
ds <- par
1832-
ds <- ds - mean(ds)
1833-
ret <- ds[index-ncat]
1834-
}
1835-
ret
1786+
as <- par[-1] * par[1]
1787+
as <- as - mean(as)
1788+
as[index]
1789+
}
1790+
fns[[i+ncat]] <- function(par, index, opar){
1791+
ds <- par - mean(par)
1792+
ds[index]
18361793
}
18371794
}
18381795
delta_index <- vector('list', ncat*2)
@@ -1846,38 +1803,21 @@ mirt2traditional <- function(x, vcov, nfact){
18461803
names(par) <- c(paste0('a', 1:ncat), paste0('c', 1:ncat))
18471804
x@est <- rep(TRUE, ncat*2)
18481805
x@SEpar <- rep(as.numeric(NA), ncat*2)
1806+
index <- as.list(c(1:ncat, 1:ncat)) ## this differs as everything changes
18491807
} else if(cls == 'nestlogit'){
18501808
if(nfact > 1L) return(x)
1851-
fns <- vector('list', ncat*2 + 4)
1852-
fns[[2]] <- function(par, index, opar){
1853-
if(index == 2L){
1854-
opar[1L:2L] <- par
1855-
ret <- -opar[2L]/opar[1L]
1856-
}
1857-
ret
1858-
}
1859-
fns[[3]] <- function(par, index, opar){
1860-
if(index == 3L)
1861-
ret <- plogis(par)
1862-
ret
1863-
}
1864-
fns[[4]] <- function(par, index, opar){
1865-
if(index == 4L)
1866-
ret <- plogis(par)
1867-
ret
1868-
}
1869-
for(i in (2L:length(par)-1L) + 4){
1809+
fns <- vector('list', (ncat-1)*2 + 4)
1810+
fns[[2]] <- function(par, index, opar) -opar[2L]/opar[1L]
1811+
fns[[3]] <- function(par, index, opar) plogis(par)
1812+
fns[[4]] <- function(par, index, opar) plogis(par)
1813+
for(i in 1:(ncat-1) + 4){
18701814
fns[[i]] <- function(par, index, opar){
1871-
if(index <= (5 + ncat-2)){
1872-
as <- par
1873-
as <- as - mean(as)
1874-
ret <- as[index-4]
1875-
} else {
1876-
ds <- par
1877-
ds <- ds - mean(ds)
1878-
ret <- ds[index-4-(ncat-1)]
1879-
}
1880-
ret
1815+
as <- par - mean(par)
1816+
as[index]
1817+
}
1818+
fns[[i+ncat-1]] <- function(par, index, opar){
1819+
ds <- par - mean(par)
1820+
ds[index]
18811821
}
18821822
}
18831823

@@ -1900,51 +1840,37 @@ mirt2traditional <- function(x, vcov, nfact){
19001840
par <- c(par1, as, ds)
19011841
x@est <- c(x@est[1:4], rep(TRUE, (ncat-1)*2))
19021842
x@SEpar <- c(x@SEpar[1], as.numeric(rep(NA, (ncat-1)*2 + 3)))
1843+
index <- delta_index
1844+
index[1:(ncat-1) + 4] <- 1:(ncat-1)
1845+
index[1:(ncat-1) + 4 + (ncat-1)] <- 1:(ncat-1)
19031846
} else if(cls == 'Luo2001'){
19041847
fns <- vector('list', nfact + ncat)
1905-
fns[[nfact+1L]] <- function(par, index, opar){
1906-
if(index == (nfact + 1L)){
1907-
opar[c(which.a, nfact + 1L)] <- par
1908-
ret <- -opar[nfact + 1L]/opar[which.a]
1909-
}
1910-
ret
1911-
}
1848+
fns[[nfact+1L]] <- function(par, index, opar) -par[2]/par[1]
19121849
for(i in 1:(ncat-1))
1913-
fns[[nfact + 1L + i]] <- function(par, index, opar){
1914-
if(index > (nfact + 1L))
1915-
ret <- exp(par)
1916-
ret
1917-
}
1850+
fns[[nfact + 1L + i]] <- function(par, index, opar) exp(par)
19181851
delta_index <- c(as.list(rep(NA, nfact)),
19191852
list(c(which.a, nfact + 1L)),
19201853
as.list((nfact + 2L):length(par)))
19211854
par[nfact + 1L] <- -par[nfact + 1L]/par[which.a]
19221855
par[(nfact + 2L):length(par)] <- exp(par[(nfact + 2L):length(par)])
19231856
names(par) <- c(a.nms, 'b', paste0('rho', 1:(ncat-1)))
1857+
index <- delta_index
19241858
}
19251859
x@par <- par
19261860
names(x@est) <- names(par)
19271861
x@parnames <- names(x@par)
19281862
if(length(vcov) == 0L || (is.na(vcov[1L,1L]) || !(cls %in%
1929-
c('dich', 'graded', 'gpcm', 'nominal', 'nestlogit')))){
1863+
c('dich', 'graded', 'gpcm', 'nominal', 'nestlogit', 'Luo2001')))){
19301864
x@SEpar <- numeric()
19311865
} else {
1932-
nms <- colnames(vcov)
1933-
splt <- strsplit(nms, "\\.")
1934-
splt <- lapply(splt, function(x) as.integer(x[-1L]))
1866+
vpick <- subset_vcov(xorg, vcov)
19351867
for(i in 1L:length(delta_index)){
19361868
if(!x@est[i]) next
19371869
if(is.na(delta_index[[i]][1L])) next
1938-
vpick <- vcov[delta_index[[i]],delta_index[[i]],drop=FALSE]
1939-
grad <- numerical_deriv(opar[delta_index[[i]]], fns[[i]], opar=opar, index=i)
1940-
parnum <- x@parnum[delta_index[[i]]]
1941-
pick <- numeric(length(grad))
1942-
for(j in 1L:length(parnum))
1943-
pick[j] <- suppressWarnings(min(which(do.call(c, lapply(splt, function(x)
1944-
any(x %in% parnum[j]))))))
1945-
grad <- grad[is.finite(pick)]
1946-
pick <- pick[is.finite(pick)]
1947-
x@SEpar[i] <- as.vector(sqrt(grad %*% vcov[pick, pick, drop=FALSE] %*% grad))
1870+
svpick <- vpick[delta_index[[i]], delta_index[[i]], drop=FALSE]
1871+
spar <- opar[delta_index[[i]]]
1872+
x@SEpar[i] <- DeltaMethod(fns[[i]], spar, svpick,
1873+
index=index[[i]], opar=opar)$se
19481874
}
19491875
if(cls == 'gpcm') x@SEpar <- x@SEpar[1L:length(x@par)]
19501876
}

‎R/utils.R‎

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -2897,8 +2897,10 @@ subset_vcov <- function(obj, vcov){
28972897
nms <- colnames(vcov)
28982898
splt <- strsplit(nms, "\\.")
28992899
splt <- lapply(splt, function(x) as.integer(x[-1L]))
2900-
pick <- sapply(splt, \(ind) ind %in% obj@parnum)
2901-
vcov[pick, pick, drop=FALSE]
2900+
pick <- sapply(splt, \(ind) any(ind %in% obj@parnum))
2901+
ret <- matrix(0, length(obj@par), length(obj@par))
2902+
ret[obj@est, obj@est] <- vcov[pick, pick, drop=FALSE]
2903+
ret
29022904
}
29032905

29042906
makeSymMat <- function(mat){

0 commit comments

Comments
 (0)