R/indices.R

In ErikOSorensen/pppindices: Calculate Standard Prices Indices For International Comparisons

#' Calculate the Törnqvist (pairwise) quantity indexes.
#' (*Not* per capita)
#'
#' @param W A ncountries x ngoods array of budget shares
#' @param W A ncountries x ngoods array of quantities
#' @return A ncountries x ncountries array of quantity indices
#'
#' @examples
#' w <- matrix( c(0.6, 0.4, 0.3, 0.4, 0.6, 0.7), 3, 2)
#' q <- matrix( c(1,2,3,4,5,6), 3, 2)
#' tornqvistq(w,q)
tornqvistq <- function(W, Q) {
  stopifnot( nrow(W)==nrow(Q), ncol(W)==ncol(Q) )
  ncountries <- nrow(W)
  ngoods <- ncol(W)
  logQt <- matrix(0, ncountries, ncountries)
  for ( i in 1:ncountries) {
    for (j in 1:ncountries) {
      for (k in 1:ngoods) {
        logQt[i,j]  <- logQt[i,j] + 0.5*(W[i,k] + W[j,k])*(log(Q[i,k]) - log(Q[j,k]))
      }
    }
  }
  exp(logQt)
}

#' Calculate the CCD quantity index.
#' The quantity index for a set of countries (*not* per capita).
#'
#' @param W A ncountries x ngoods array budget shares
#' @param Q A ncountries x ngoods array of quantities
#' @return A ncountries long vector of the CCD quantity CCD[j], relative to minimum income
#'
#' @examples
#' w <- matrix( c(0.6, 0.4, 0.3, 0.4, 0.6, 0.7), 3, 2)
#' q <- matrix( c(1,2,3,4,5,6), 3, 2)
#' CCDq(w,q)

CCDq <- function(W, Q) {
  stopifnot( nrow(W)==nrow(Q), ncol(W)==ncol(Q) )
  ncountries <- nrow(W)
  ngoods <- ncol(W)
  logTq <- log( tornqvistq(W, Q))
  logCCD <- rowMeans(logTq)
  exp(logCCD - min(logCCD))
}

#' Calculate the Fisher Ideal quantity index
#' (*not* per capita).
#'
#' @param P A ncountries x ngoods array of prices.
#' @param Q A ncountries x ngoods array of quantities
#'
#' @examples
#' p <- matrix( c(0.6, 0.4, 0.3, 0.4, 0.6, 0.7), 3, 2)
#' q <- matrix( c(1,2,3,4,5,6), 3, 2)
#' fisheridealq(w,q)
fisheridealq <- function(P, Q) {
  stopifnot( nrow(P)==nrow(Q), ncol(P)==ncol(Q) )
  ncountries = ncol(P)
  ngoods = ncol(P)
  logQF <- matrix(0, ncol=ncountries, nrow=ncountries)
  for (i in 1:ncountries) {
    for (j in 1:ncountries) {
      0.5 * log( P[,j] %*% Q[,i] / P[,j] %*% Q[,j] ) +
      0.5 * log( P[,i] %*% Q[,i] / P[,i] %*% Q[,j] )
    }
  }
  exp(logQF)
}