Create functions

# ICC 4p
       
        irt_4p <- function(theta,a,b,c,d){c + (d-c)*(exp(a*(theta-b))/(1+(exp(a*(theta-b)))))}

# Information function 4p
        info_4p <- function(theta,a,b,c,d){
                a_sqrd <- a**2
                p_theta <- irt4p(theta, a,b,c,d)
    
                i <- ((a_sqrd*((p_theta-c)**2)*((d - p_theta)**2)) 
                / 
                (((d-c)**2)*p_theta*(1-p_theta)))
                
          return(i)
        }

# Optimumm weight 4p
# See Magis (2013) p. 306 e Lord (1980) p.74
# Note that in MIRT notation (c=g e u=d). We are using traditional a, b, c and d notation.
       
        optm_weight_4p <- function(theta,a,b,c,d){
                
                p_theta <- irt_4p(theta, a,b,c,d)
                
                i <- (
                        (a / (d - c))*(p_theta - c)*(d - p_theta)
                        / 
                        (p_theta*(1-p_theta))
                   )
          return(i)
        }

Example 1 (3-parameter model): item parameter’s, response pattern’s and data

# Item parameters
        bs <- c(-2.4, -2.0, -2.0, -1.0, 0, 1.2, 2.0, 2.0, 2.2, 2.4)
        as <- c(0.8, 0.8, 0.8, 0.8, 0.8, 1.2, 1.2, 1.2, 1.2, 1.2)
        
        set.seed(44)
        cs <- rbeta(n = 10, shape1 = 5, shape2 = 17)
        round(cs, digits = 2)
##  [1] 0.30 0.23 0.08 0.21 0.12 0.11 0.34 0.21 0.09 0.15
        ds <- rep(1, 10)

# Thetas         
        theta <- seq(from=-4, to=4, by = 1)
        
        data <- as.data.frame(theta)
        
# Data 
        data <-  cbind(data, 
                lapply(
                1:10, function(x)
                {
                  round(optm_weight_4p(theta, as[x], bs[x], cs[x], ds[x]),3)
                }
          )
        )
        
        names(data) [2:11] <- paste("i", 1:10, sep = "_")
      
# Response Patterns
        r_ptrn1 <-c(1,1,1,1,0,0,0,0,0,0)
        r_ptrn2 <-c(0,0,0,0,0,0,1,1,1,1)
        library(catR)
## Warning: package 'catR' was built under R version 3.3.2
        it_param <- as.matrix( cbind(as, bs, cs, ds) )
        row_names <- (names(data)[2:11])
        col_names <- c("a", "b", "c", "d")
        dimnames(it_param) <- list(row_names, col_names)

# Latent scores dor each pattenr       
        thetaEst(it = it_param, x = r_ptrn1, method = "ML")
## [1] -0.5469524
        thetaEst(it = it_param, x = r_ptrn2, method = "ML")
## [1] -4
# matrix of optimum weights        
        library(tidyr)
## Warning: package 'tidyr' was built under R version 3.3.2
        library(ggplot2)
## Warning: package 'ggplot2' was built under R version 3.3.2
        data %>% 
          gather(2:11, key = "items", value = "weights") %>%      
          ggplot(aes(y = theta, x = factor(items, levels = row_names))) +
           geom_tile(aes(fill = weights, alpha = .4)) + 
           geom_text(aes(label = round(weights, 1)), size=4) +
            theme_bw() +
           scale_fill_gradient(low = "white", high = "orange") +
           scale_y_continuous(trans = "reverse", breaks = theta) +
           theme(axis.title.x=element_blank())

Example 2 (4-parameter model): item parameter’s, response pattern’s and data

# Item parameters
        bs <- c(-2.4, -2.0, -2.0, -1.0, 0, 1.2, 2.0, 2.0, 2.2, 2.4)
        as <- c(0.8, 0.8, 0.8, 0.8, 0.8, 1.2, 1.2, 1.2, 1.2, 1.2)        
        cs <- rep(0, 10)
        ds <- c(.6, .6, .6, 1, 1, 1, 1, 1, .6, .6)

# Theta        
        theta <- seq(from=-4, to=4, by = 1)

# Data        
        data2 <- as.data.frame(theta)
        
        data2 <-  cbind(data2, 
                lapply(
                1:10, function(x)
                {
                  round(optm_weight_4p(theta, as[x], bs[x], cs[x], ds[x]),3)
                }
          )
        )
        names(data2) [2:11] <- paste("i", 1:10, sep = "_")
       
# Response patterns
        r_ptrn1 <-c(0,0,1,1,1,1,1,1,1,0)
        r_ptrn2 <-c(1,1,1,1,1,1,1,1,1,0)

        
# Theta estimated for each response pattenr     
        library(catR)
        
        it_param <- as.matrix( cbind(as, bs, cs, ds) )
        row_names <- (names(data2)[2:11])
        col_names <- c("a", "b", "c", "d")
        dimnames(it_param) <- list(row_names, col_names)
       
        thetaEst(it = it_param, x = r_ptrn1, method = "ML")
## [1] 4
        thetaEst(it = it_param, x = r_ptrn2, method = "ML")
## [1] 4
# Matrix of optimum weights

        library(tidyr)
        data2 %>% 
          gather(2:11, key = "items", value = "weights") %>%      
          ggplot(aes(y = theta, x = factor(items, levels = row_names))) +
           geom_tile(aes(fill = weights, alpha = .4)) + 
           geom_text(aes(label = round(weights, 1))) +
            theme_bw() +
           scale_fill_gradient(low = "white", high = "orange") +
           scale_y_continuous(trans = "reverse", breaks = theta) +
           theme(axis.title.x=element_blank())