[Eventstudies-commits] r73 - in pkg: R data man

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Sat Jul 13 12:30:57 CEST 2013


Author: vikram
Date: 2013-07-13 12:30:57 +0200 (Sat, 13 Jul 2013)
New Revision: 73

Added:
   pkg/R/AMM.R
   pkg/R/firmExposures.R
   pkg/R/makeX.R
   pkg/R/manyfirmsAMM.R
   pkg/R/onefirmAMM.R
   pkg/data/ammData.rda
   pkg/data/firmExposuresData.rda
   pkg/man/AMM.Rd
   pkg/man/ammData.Rd
   pkg/man/firmExposures.Rd
   pkg/man/firmExposuresData.Rd
   pkg/man/makeX.Rd
   pkg/man/manyfirmsAMM.Rd
   pkg/man/onefirmAMM.Rd
Log:
Added AMM functionality function and related documentation files and codebase; Work in progress

Added: pkg/R/AMM.R
===================================================================
--- pkg/R/AMM.R	                        (rev 0)
+++ pkg/R/AMM.R	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,470 @@
+##########################
+# Generalised AMM function
+##########################
+AMM <- function(amm.type=c("onefirm","manyfirms","firmExposures"), ...){
+  #-----------------------------
+  # Getting RHS variables in AMM
+  #-----------------------------
+  print("Preparing explanatory variables for computation of Augmented Market Models")
+  print("If any argument is missing then please refer to the documentation for detailed explanation.")
+  # Checking arguments  
+  nlags.check <- testObject(nlags)
+  if(nlags.check==FALSE){ nlags <- NA}
+  verbose.check <- testObject(verbose)
+  if(verbose.check==FALSE){ verbose <- FALSE}
+  dates.check <- testObject(dates)
+  if(dates.check==FALSE){ dates <- NULL}
+  
+  rM1.check <- testObject(rM1)
+  if(rM1.check==FALSE){ stop("Input rM1 (stock market index) is missing")}
+  others.check <- testObject(others)
+  if(others.check==FALSE){ stop("Input 'others' (time series of other regressor or interest) is missing")}
+  rM1.check <- testObject(rM1)
+  if(rM1.check==FALSE){ stop("Input rM1 (stock market index) is missing")}
+  rM1purge.check <- testObject(rM1purge)
+  if(rM1purge.check==FALSE){ stop("Input rM1purge is missing")}
+  switch.to.innov.check <- testObject(switch.to.innov)
+  if(switch.to.innov.check==FALSE){ stop("Input switch.to.innov is missing")}
+
+  if(amm.type!="manyfirms"){
+    X <- makeX(rM1, others, switch.to.innov,
+               rM1purge, nlags, dates, verbose)
+  }else{
+    regressors <- makeX(rM1, others, switch.to.innov, rM1purge, nlags,
+                        dates, verbose)
+  }
+  
+  #---------
+  # One firm
+  #---------
+  if(amm.type=="onefirm"){
+    # Checking if arguments are provided
+    rj.check <- testObject(rj)
+    if(rj.check==FALSE){ stop("Input rj (firm data) is missing")}
+ 
+    result <- onefirmAMM(rj, X, nlags, verbose, dates)
+  }
+    
+  #-----------
+  # Many firms
+  #-----------
+  if(amm.type=="manyfirms"){
+    # Checking if arguments are provided
+    regressors.check <- testObject(regressors)
+    if(regressors.check==FALSE){ stop("Input regressors is missing. Refer documentation.")}
+    regressand.check <- testObject(regressand)
+    if(regressand.check==FALSE){ stop("Input regressand is missing. Refer documentation.")}
+    dates.check <- testObject(dates)
+    if(dates.check==FALSE){ dates <- NULL}
+    periodnames.check <- testObject(periodnames)
+    if(periodnames.check==FALSE){ periodnames <- NULL}
+    
+    result <- manyfirmsAMM(regressand,regressors,
+                        lags=nlags,dates, periodnames,verbose)
+  }
+
+  #---------------
+  # Firm exposures
+  #---------------
+  if(amm.type=="firmExposures"){
+    rj.check <- testObject(rj)
+    if(rj.check==FALSE){ stop("Input rj (firm data) is missing")}
+    
+    result <- firmExposures(rj, X, nlags, verbose)
+  }
+
+  return(result)
+}
+
+#######################
+# AMM for one firm
+#######################
+onefirmAMM <- function(rj,X,nlags=NA,verbose=FALSE,dates=NULL){
+  exposures <- data.frame(matrix(NA,ncol=ncol(X),nrow=(length(dates)-1)))
+  colnames(exposures) <- colnames(X)
+  sds <- exposures
+  periodnames <- NULL
+  
+  if(is.null(dates)){
+   res <- firmExposures(rj,X,verbose=verbose,nlags=nlags)
+   exposures <- res$exposure
+   sds <- res$s.exposure
+ }else{
+   for(i in 1:(length(dates)-1)){
+     tmp <- window(rj,start=dates[i],end=dates[i+1])
+     rhs <- window(X,start=dates[i],end=dates[i+1])
+     res <- firmExposures(rj=tmp,
+                          X=rhs,
+                          verbose=verbose,
+                          nlags=nlags)
+     exposures[i,] <- res$exposure
+     periodnames <- c(periodnames,paste(dates[i],dates[i+1],sep=" TO "))
+     sds[i,] <- res$s.exposure
+   }
+   rownames(exposures) <- rownames(sds) <- periodnames
+ }
+  rval <- list(exposures=exposures,sds=sds)
+  return(rval)
+}
+
+########################
+# Many firms AMM
+########################
+manyfirmsAMM <-function(regressand,regressors,
+                        lags,dates=NULL, periodnames=NULL,verbose=FALSE){
+  if(is.null(dates)){
+    dates=c(start(regressors),end(regressors))
+    periodnames="Full"
+  }
+  nperiods <- length(periodnames)
+  if(length(dates) != (nperiods+1)){
+    cat("Mistake in length of dates versus length of periods.\n")
+    return(NULL)
+  }
+
+  nfirms <- ncol(regressand)
+
+  # Let's get "exposure' and 'sds'. Setting up structures:-
+  exposures <- matrix(NA,nrow=nfirms,ncol=nperiods*ncol(regressors))
+  rownames(exposures) <- colnames(regressand)
+  tmp <- NULL
+  for(i in 1:length(periodnames)){
+    for(j in 1:ncol(regressors)){
+      tmp <-  c(tmp, paste(colnames(regressors)[j],
+                           periodnames[i],sep="."))
+    }
+  }
+  colnames(exposures) <- tmp
+  sds <- exposures
+  colnames(sds) <- paste("sd",colnames(exposures),sep=".")
+
+  # Setup a list structure for an OLS that failed
+  empty <- list(exposures=rep(NA,ncol(regressors)),
+                s.exposures=rep(NA,ncol(regressors)))
+
+  for(i in 1:ncol(regressand)){
+	cat("Doing",colnames(regressand)[i])
+    if (verbose) {cat ("Doing", colnames(regressand)[i])}
+    rj <- regressand[,i]
+    dataset <- cbind(rj, regressors)   # This is the full time-series
+    this.exp <- this.sds <- NULL
+    for(j in 1:nperiods){              # now we chop it up 
+      t1 <- dates[j]
+      t2 <- dates[j+1]
+     # if (j != nperiods) {t2 <- t2-1} # I don't know why I wrote this now :-( But this created problems.
+      this <- window(dataset,start=t1, end=t2)
+      fe <- firmExposures(this[,1],this[,-1],nlags=lags,verbose)
+      if(is.null(fe)) {fe <- empty}
+      this.exp <- c(this.exp, fe$exposures)
+      this.sds <- c(this.sds, fe$s.exposures)
+    }
+    exposures[i,] <- this.exp
+    sds[i,] <- this.sds
+  }
+  list(exposures=exposures, sds=sds, sig=exposures/sds)
+}
+
+###############################################
+# Estimating one firm's exposure in one period.
+###############################################
+firmExposures <- function(rj, X, nlags=NA, verbose=FALSE) {
+  do.ols <- function(nlags) {
+    tmp <- cbind(rj, X[,1]) # Assume 1st is stock index, and no lags are required there.
+    labels <- c("rj","rM1")
+    if (NCOL(X) > 1) {
+      for (i in 2:NCOL(X)) {
+        for (j in 0:nlags) {
+          tmp <- cbind(tmp, lag(X[,i], -j))
+          labels <- c(labels, paste(colnames(X)[i], j, sep="."))
+        }
+      }
+    }
+    tmp <- na.omit(tmp)
+    if (nrow(tmp) < 30) {             # refuse to do the work.
+      return(NULL)                    # returns out of do.ols() only
+    }
+
+    colnames(tmp) <- labels          # So the OLS results will look nice
+    lm(rj ~ ., data=as.data.frame(tmp))
+  }
+
+  if (is.na(nlags)) {
+    if (verbose) {cat("Trying to find the best lag structure...\n")}
+    bestlag <- 0
+    bestm <- NULL
+    bestAIC <- Inf
+    for (trylag in 0:min(10,log10(length(rj)))) {
+      thism <- do.ols(trylag)
+      thisAIC <- AIC(thism, k=log(length(thism$fitted.values)))
+      if (verbose) {cat(trylag, " lags, SBC = ", thisAIC, "\n")}
+      if (thisAIC < bestAIC) {
+        bestlag <- trylag
+        bestAIC <- thisAIC
+        bestm <- thism
+      }
+    }
+    nlags <- bestlag
+    m <- bestm
+  } else {
+    m <- do.ols(nlags)
+    if (is.null(m)) {return(NULL)}
+  }
+  # In either event, you endup holding an "m" here.
+  if (verbose) {cat("\n\nThe OLS:\n"); print(summary(m))}
+
+  # Compute a series of exposure measures, and their standard errors.
+  beta <- m$coefficients
+  Sigma <- vcovHAC(m)
+  # First the rM1
+  exposures <- beta[2]                  # no lags for rM1
+  s.exposures <- sqrt(Sigma[2,2])
+  # From here on, there's a block of 1+nlags coeffs for each
+  # of the non-rM1 regressors.
+  if (NCOL(X) > 1) {
+    for (i in 2:NCOL(X)) {
+      n.block1 <- 2 + ((i-2)*(1+nlags)) # Just 2 for the 1st case.
+      n.block2 <- length(beta) - n.block1 - (1 + nlags)
+      w <- c(rep(0, n.block1), rep(1, 1+nlags), rep(0, n.block2))
+      exposures <- c(exposures, w %*% beta)
+      s.exposures <- c(s.exposures, sqrt(w %*% Sigma %*% w))
+    }
+  }
+  names(exposures) <- names(s.exposures) <- colnames(X)
+  results <- list(exposures=exposures,
+                  s.exposures=s.exposures, nlags=nlags,
+                  lm.res=m)
+  class(results) <- "amm"
+  results
+}
+
+###########################
+# Maintaing NAs in AR model
+###########################
+ARinnovations <- function(x) {
+  stopifnot(NCOL(x) == 1)
+  dt <- NULL
+  if (class(x) == "zoo") {
+    dt <- index(x)
+    x <- as.numeric(x)
+  }
+  non.na.locations <- !is.na(x)
+  x <- x[non.na.locations]
+  m <- ar(x)
+  e <- m$resid
+  # Relocate with the old NA structure
+  result <- rep(NA, length(non.na.locations))
+  result[non.na.locations] <- e
+  if (!is.null(dt)) {result <- zoo(result, order.by=dt)}
+  list(result=result, m=m)
+}
+
+# ---------------------------------------------------------------------------
+# The workhorse called by makeX to return a nice matrix of RHS
+# variables to be used in an analysis. 
+do.one.piece <- function(rM1, others, switch.to.innov, rM1purge, nlags, verbose=FALSE) {
+  cat("   do.one.piece invocation\n")
+  thedates <- index(rM1)
+  if (verbose) {
+    cat("   Doing work for period from ",
+        as.character(head(thedates,1)), " to ",
+        as.character(tail(thedates,1)), "\n")
+  }
+  otherlags <- NULL
+  for (i in 1:NCOL(others)) {
+    if (switch.to.innov[i]) {
+      a <- ARinnovations(others[,i])
+      innovated <- a$result
+      if (verbose) {
+        cat("   AR model for ", colnames(others, do.NULL=FALSE)[i], "\n")
+        print(a$m)
+      }
+      otherlags <- c(otherlags, a$m$order)
+    } else {
+      innovated <- others[,i]
+      otherlags <- c(otherlags, 0)
+    }
+    if (i > 1) {
+      innov <- cbind(innov, innovated)
+    } else {
+      innov <- innovated
+    }
+  }
+  if (NCOL(innov) > 1) {colnames(innov) <- colnames(others)}
+  rM1.purged <- rM1
+  if (rM1purge) {
+    firstpass <- TRUE
+    for (i in 1:NCOL(innov)) {
+      for (j in 0:nlags) {
+        if (firstpass) {
+          z <- lag(innov[,i],-j)
+          labels <- paste(colnames(innov,do.NULL=FALSE)[i],j,sep=".")
+          firstpass <- FALSE
+        } else {
+          z <- cbind(z, lag(innov[,i], -j))
+          labels <- c(labels, paste(colnames(innov,do.NULL=FALSE)[i],j,sep="."))
+        }
+      }
+    }
+    if (NCOL(z) > 1) {colnames(z) <- labels}
+    m <- lm(rM1 ~ ., as.data.frame(cbind(rM1, z)))
+    if (verbose) {
+      cat("   Model explaining rM1:\n")
+      print(summary(m))
+    }
+    how.many.NAs <- nlags + max(otherlags)
+    rM1.purged <- zoo(c(rep(NA,how.many.NAs),m$residuals),
+                      order.by=thedates)
+  }
+                                        #    if (verbose) {cat("   Finished do.one.piece()\n")}
+  list(rM1.purged=rM1.purged, innov=innov)
+}                              
+
+# A function that calls do.one.piece, and works through several
+# different periods to provide the right RHS matrix. 
+makeX <- function(rM1, others,
+                  switch.to.innov=rep(TRUE, NCOL(others)),
+                  rM1purge=TRUE,
+                  nlags=5,
+                  dates=NULL,
+                  verbose=TRUE) {
+  if (verbose) {cat("0. Checking args\n")}
+  stopifnot(all.equal(index(rM1), index(others)),
+            length(switch.to.innov)==NCOL(others))
+  if (!is.null(dates)) {
+    stopifnot(class(dates) == "Date")
+  }
+  if (verbose) {cat("1. Checking dates.\n")}
+  if (is.null(dates)) {
+    dates <- c(start(rM1),end(rM1))
+  }
+  if(head(dates,1)!=head(index(rM1),1)){
+    stop("Start date provided and the start date of the dataset do not match \n")
+  }
+  if(tail(dates,1)!=tail(index(rM1),1)){
+    stop("End date provided and the end date of the dataset do not match \n")
+  }
+  if (verbose) {cat("2. Run through all the pieces --\n")}
+  for (i in 1:(length(dates)-1)) {
+    t1 <- dates[i]
+    t2 <- dates[i+1]
+    if (i != (length(dates)-1)) {t2 <- t2 -1}
+    if (verbose) {
+      cat("   Focusing down from date = ", as.character(t1), " to ", as.character(t2), "\n")
+    }
+    tmp.rM1 <- window(rM1, start=t1, end=t2)
+    tmp.others <- window(others, start=t1, end=t2)
+    a <- do.one.piece(tmp.rM1, tmp.others, switch.to.innov, rM1purge, nlags, verbose)
+    if (i > 1) {
+      res.rM1 <- c(res.rM1, a$rM1.purged)
+      res.innov <- rbind(res.innov, a$innov)
+    } else {
+      res.rM1 <- a$rM1.purged
+      res.innov <- a$innov
+    }
+  }
+  if (verbose) {cat("2. Make a clean X and send it back --\n")}
+  X <- cbind(res.rM1, res.innov)
+  if (NCOL(res.innov) == 1) {colnames(X) <- c("rM1","z")}
+  else {colnames(X) <- c("rM1", colnames(res.innov))}
+  X
+}
+
+#############################################################################
+# Exposure is a vector of exposures
+# var.exposures is a vector of the s.error of each exposure
+# X is a matrix of firm characteristics
+# Note: They all should have the same number of rows. 
+lm.for.amm.exposure <- function(exposures, var.exposures, X) {
+  m.ols <- lm(exposures ~ -1 + X)  
+# Likelihood function for OLS-by-MLE where each obs has a known standard error --
+  mine.lf <- function(theta, y, X, y.se2) {
+    if (theta[1] <= 0) return(NA)         # Refuse to look at -ve sigma.
+    -sum(dnorm(y, mean= X %*% theta[-1], sd=sqrt(theta[1]+y.se2), log=TRUE))
+  }  
+ ## Use OLS as starting values for the MLE --
+  seeds <- c(summary(m.ols)$sigma, coef(m.ols))
+  p <- optim(seeds, method="L-BFGS-B", fn=mine.lf, hessian=TRUE,
+             lower=c(1e-2, -Inf, -Inf, -Inf), upper=rep(Inf, 4),
+             y=exposures, X=X, y.se2=var.exposures)
+  inverted <- solve(p$hessian)
+ # Now get to generation of a nicely organised matrix collecting up the results.
+ # First setup a clean set of OLS results --
+  results <- summary(m.ols)$coefficients[,c(1,3)]
+  rownames(results) <- substr(rownames(results), 2, 1000)
+  colnames(results) <- c("OLS estimates", "OLS t stats")
+  results <- rbind(results, c(summary(m.ols)$sigma, "NA"))
+  rownames(results)[nrow(results)] <- "Sigma"
+# Now augment it on the right with the MLE --
+  tmp <- p$par/sqrt(diag(inverted))
+  results <- cbind(results, cbind(c(p$par[-1], p$par[1]),
+                                  c(tmp[-1], tmp[1])))
+  colnames(results) <- c(colnames(results)[1:2], "MLE estimates", "MLE t stats")
+  results
+}
+
+
+####################################################
+
+print.amm <-
+function(amm, verbose=FALSE) {
+  if (verbose) {print(summary(amm$lm.res))}
+  cat("Summary statistics of exposure:\n")
+  sstats <- cbind(amm$exposure, amm$s.exposure,
+                  amm$exposure/amm$s.exposure)
+  colnames(sstats) <- c("Exposure", "Std.Err", "t statistic")  
+  rownames(sstats) <- names(amm$exposures)
+  print(sstats)
+  return(0)
+}
+
+################################################################################
+# You got to write a plot function that does cool stuff with the AMM results! 
+################################################################################
+
+simulated.mean <- function (beta, stdev, fun=NULL, no.of.draws){
+    na1 <- is.na(beta)
+    na2 <- is.na(stdev)
+   if (!identical(na1, na2)) {
+        cat("Panic - two is.na() are not identical.\n")
+       return(-1)
+   }
+    b <- beta[!is.na(beta)]
+   s <- stdev[!is.na(stdev)]
+    N <- length(b)
+    draws <- NULL
+    for (i in 1:no.of.draws) {
+     if(!is.null(fun)){
+        draws <- c(draws, mean(fun((rnorm(N) * s) + b)))
+      }else{
+       draws <- c(draws, mean(rnorm(N)*s+b))
+     }
+  }
+    draws
+}
+
+#########################################
+kernel.plots <- function(draws, logscale=NULL) {
+    nplots <- ncol(draws)
+      hilo <- range(draws)
+      par(mfrow=c(nplots,1), mai=c(.4,.8,.2,.2))
+      for (i in 1:nplots) {
+        if(!is.null(logscale)){
+            plot(density(draws[,i]), main=colnames(draws)[i], xlab="",
+                          col="blue", xlim=hilo, log=logscale,lwd=2)
+          }else{
+	    plot(density(draws[,i]), main=colnames(draws)[i], xlab="",
+                          col="blue", xlim=hilo, lwd=2)
+          }
+	}
+  }
+
+#############################################################################################
+
+
+
+###########################
+# Checking if object exists
+###########################
+testObject <- function(object){
+   exists(as.character(substitute(object)))
+}

Added: pkg/R/firmExposures.R
===================================================================
--- pkg/R/firmExposures.R	                        (rev 0)
+++ pkg/R/firmExposures.R	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,73 @@
+###############################################
+# Estimating one firm's exposure in one period.
+###############################################
+firmExposures <- function(rj, X, nlags=NA, verbose=FALSE) {
+  do.ols <- function(nlags) {
+    tmp <- cbind(rj, X[,1]) # Assume 1st is stock index, and no lags are required there.
+    labels <- c("rj","rM1")
+    if (NCOL(X) > 1) {
+      for (i in 2:NCOL(X)) {
+        for (j in 0:nlags) {
+          tmp <- cbind(tmp, lag(X[,i], -j))
+          labels <- c(labels, paste(colnames(X)[i], j, sep="."))
+        }
+      }
+    }
+    tmp <- na.omit(tmp)
+    if (nrow(tmp) < 30) {             # refuse to do the work.
+      return(NULL)                    # returns out of do.ols() only
+    }
+
+    colnames(tmp) <- labels          # So the OLS results will look nice
+    lm(rj ~ ., data=as.data.frame(tmp))
+  }
+
+  if (is.na(nlags)) {
+    if (verbose) {cat("Trying to find the best lag structure...\n")}
+    bestlag <- 0
+    bestm <- NULL
+    bestAIC <- Inf
+    for (trylag in 0:min(10,log10(length(rj)))) {
+      thism <- do.ols(trylag)
+      thisAIC <- AIC(thism, k=log(length(thism$fitted.values)))
+      if (verbose) {cat(trylag, " lags, SBC = ", thisAIC, "\n")}
+      if (thisAIC < bestAIC) {
+        bestlag <- trylag
+        bestAIC <- thisAIC
+        bestm <- thism
+      }
+    }
+    nlags <- bestlag
+    m <- bestm
+  } else {
+    m <- do.ols(nlags)
+    if (is.null(m)) {return(NULL)}
+  }
+  # In either event, you endup holding an "m" here.
+  if (verbose) {cat("\n\nThe OLS:\n"); print(summary(m))}
+
+  # Compute a series of exposure measures, and their standard errors.
+  beta <- m$coefficients
+  Sigma <- vcovHAC(m)
+  # First the rM1
+  exposures <- beta[2]                  # no lags for rM1
+  s.exposures <- sqrt(Sigma[2,2])
+  # From here on, there's a block of 1+nlags coeffs for each
+  # of the non-rM1 regressors.
+  if (NCOL(X) > 1) {
+    for (i in 2:NCOL(X)) {
+      n.block1 <- 2 + ((i-2)*(1+nlags)) # Just 2 for the 1st case.
+      n.block2 <- length(beta) - n.block1 - (1 + nlags)
+      w <- c(rep(0, n.block1), rep(1, 1+nlags), rep(0, n.block2))
+      exposures <- c(exposures, w %*% beta)
+      s.exposures <- c(s.exposures, sqrt(w %*% Sigma %*% w))
+    }
+  }
+  names(exposures) <- names(s.exposures) <- colnames(X)
+  results <- list(exposures=exposures,
+                  s.exposures=s.exposures, nlags=nlags,
+                  lm.res=m)
+  class(results) <- "amm"
+  results
+}
+

Added: pkg/R/makeX.R
===================================================================
--- pkg/R/makeX.R	                        (rev 0)
+++ pkg/R/makeX.R	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,54 @@
+
+###############################################################
+# A function that calls do.one.piece, and works through several
+# different periods to provide the right RHS matrix.
+###############################################################
+source('AMM.R')
+makeX <- function(rM1, others,
+                  switch.to.innov=rep(TRUE, NCOL(others)),
+                  rM1purge=TRUE,
+                  nlags=5,
+                  dates=NULL,
+                  verbose=TRUE) {
+  if (verbose) {cat("0. Checking args\n")}
+  stopifnot(all.equal(index(rM1), index(others)),
+            length(switch.to.innov)==NCOL(others))
+  if (!is.null(dates)) {
+    stopifnot(class(dates) == "Date")
+  }
+  if (verbose) {cat("1. Checking dates.\n")}
+  if (is.null(dates)) {
+    dates <- c(start(rM1),end(rM1))
+  }
+  if(head(dates,1)!=head(index(rM1),1)){
+    stop("Start date provided and the start date of the dataset do not match \n")
+  }
+  if(tail(dates,1)!=tail(index(rM1),1)){
+    stop("End date provided and the end date of the dataset do not match \n")
+  }
+  if (verbose) {cat("2. Run through all the pieces --\n")}
+  for (i in 1:(length(dates)-1)) {
+    t1 <- dates[i]
+    t2 <- dates[i+1]
+    if (i != (length(dates)-1)) {t2 <- t2 -1}
+    if (verbose) {
+      cat("   Focusing down from date = ", as.character(t1), " to ", as.character(t2), "\n")
+    }
+    tmp.rM1 <- window(rM1, start=t1, end=t2)
+    tmp.others <- window(others, start=t1, end=t2)
+    a <- do.one.piece(tmp.rM1, tmp.others, switch.to.innov, rM1purge, nlags, verbose)
+    if (i > 1) {
+      res.rM1 <- c(res.rM1, a$rM1.purged)
+      res.innov <- rbind(res.innov, a$innov)
+    } else {
+      res.rM1 <- a$rM1.purged
+      res.innov <- a$innov
+    }
+  }
+  if (verbose) {cat("2. Make a clean X and send it back --\n")}
+  X <- cbind(res.rM1, res.innov)
+  if (NCOL(res.innov) == 1) {colnames(X) <- c("rM1","z")}
+  else {colnames(X) <- c("rM1", colnames(res.innov))}
+  X
+}
+

Added: pkg/R/manyfirmsAMM.R
===================================================================
--- pkg/R/manyfirmsAMM.R	                        (rev 0)
+++ pkg/R/manyfirmsAMM.R	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,58 @@
+
+########################
+# Many firms AMM
+########################
+source('AMM.R')
+manyfirmsAMM <-function(regressand,regressors,
+                        lags,dates=NULL, periodnames=NULL,verbose=FALSE){
+  if(is.null(dates)){
+    dates=c(start(regressors),end(regressors))
+    periodnames="Full"
+  }
+  nperiods <- length(periodnames)
+  if(length(dates) != (nperiods+1)){
+    cat("Mistake in length of dates versus length of periods.\n")
+    return(NULL)
+  }
+  nfirms <- ncol(regressand)
+
+  # Let's get "exposure' and 'sds'. Setting up structures:-
+
+  exposures <- matrix(NA,nrow=nfirms,ncol=nperiods*ncol(regressors))
+  rownames(exposures) <- colnames(regressand)
+  tmp <- NULL
+  for(i in 1:length(periodnames)){
+    for(j in 1:ncol(regressors)){
+      tmp <-  c(tmp, paste(colnames(regressors)[j],
+                           periodnames[i],sep="."))
+    }
+  }
+  colnames(exposures) <- tmp
+  sds <- exposures
+  colnames(sds) <- paste("sd",colnames(exposures),sep=".")
+
+  # Setup a list structure for an OLS that failed
+  empty <- list(exposures=rep(NA,ncol(regressors)),
+                s.exposures=rep(NA,ncol(regressors)))
+
+  for(i in 1:ncol(regressand)){
+	cat("Doing",colnames(regressand)[i])
+    if (verbose) {cat ("Doing", colnames(regressand)[i])}
+    rj <- regressand[,i]
+    dataset <- cbind(rj, regressors)   # This is the full time-series
+    this.exp <- this.sds <- NULL
+    for(j in 1:nperiods){              # now we chop it up 
+      t1 <- dates[j]
+      t2 <- dates[j+1]
+     # if (j != nperiods) {t2 <- t2-1} # I don't know why I wrote this now :-( But this created problems.
+      this <- window(dataset,start=t1, end=t2)
+      fe <- firmExposures(this[,1],this[,-1],nlags=lags,verbose)
+      if(is.null(fe)) {fe <- empty}
+      this.exp <- c(this.exp, fe$exposures)
+      this.sds <- c(this.sds, fe$s.exposures)
+    }
+    exposures[i,] <- this.exp
+    sds[i,] <- this.sds
+  }
+  list(exposures=exposures, sds=sds, sig=exposures/sds)
+}

Added: pkg/R/onefirmAMM.R
===================================================================
--- pkg/R/onefirmAMM.R	                        (rev 0)
+++ pkg/R/onefirmAMM.R	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,31 @@
+#######################
+# AMM for one firm
+#######################
+source('AMM.R')
+onefirmAMM <- function(rj,X,nlags=NA,verbose=FALSE,dates=NULL){
+  exposures <- data.frame(matrix(NA,ncol=ncol(X),nrow=(length(dates)-1)))
+  colnames(exposures) <- colnames(X)
+  sds <- exposures
+  periodnames <- NULL
+  
+  if(is.null(dates)){
+   res <- firmExposures(rj,X=rhs,verbose=verbose,nlags=nlags)
+   exposures <- res$exposure
+   sds <- res$s.exposure
+ }else{
+   for(i in 1:(length(dates)-1)){
+     tmp <- window(rj,start=dates[i],end=dates[i+1])
+     rhs <- window(X,start=dates[i],end=dates[i+1])
+     res <- firmExposures(rj=tmp,
+                          X=rhs,
+                          verbose=verbose,
+                          nlags=nlags)
+     exposures[i,] <- res$exposure
+     periodnames <- c(periodnames,paste(dates[i],dates[i+1],sep=" TO "))
+     sds[i,] <- res$s.exposure
+   }
+   rownames(exposures) <- rownames(sds) <- periodnames
+ }
+  rval <- list(exposures=exposures,sds=sds)
+  return(rval)
+}

Added: pkg/data/ammData.rda
===================================================================
(Binary files differ)


Property changes on: pkg/data/ammData.rda
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/data/firmExposuresData.rda
===================================================================
(Binary files differ)


Property changes on: pkg/data/firmExposuresData.rda
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/man/AMM.Rd
===================================================================
--- pkg/man/AMM.Rd	                        (rev 0)
+++ pkg/man/AMM.Rd	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,67 @@
+\name{AMM}
+\alias{AMM}
+
+\title{A function that estimates exposure for a single firm over multiple periods}
+
+\description{This function does three tasks: One firm AMM model, many firms AMM model and firm exposures. One should refer to the individual functions manual for arguments. 
+}
+
+\usage{
+AMM(amm.type=c("one.firm","manyfirms","firmExposures"), ...)
+}
+
+\arguments{
+  \item{amm.type}{Type of AMM to be performed: On one firm or many firms or just firm exposure}
+
+}
+\value{ The function returns the exposures, HAC adjusted standard
+  errors, the number of lags used, and the residuals from the fitted
+  model.
+
+\item{exposures}{This contains the exposure estimates for the firm
+j}
+
+\item{s.exposures}{This contains the HAC adjusted standard error
+of the exposures estimated for the firm rj.}
+
+\item{nlags}{Specifies the lag length provided by the user}
+
+\item{lm.res}{The model estimates used in the analysis are stored here.}
+}
+
+\author{Vikram Bahure}
+
+\seealso{ \code{\link{firmExposures}},
+\code{\link{manyfirmsAMM}}
+\code{\link{onefirmAMM}}}
+
+\examples{ 
+# Create RHS before running AMM()
+data("ammData",package="amm")
+NIFTY_INDEX <- y3c3$NIFTY_INDEX
+INRUSD <- y3c3$INRUSD
+Company_A <- y3c3$Company_A
+Company_B <- y3c3$Company_B
+Company_C <- y3c3$Company_C
+regressand <- cbind(Company_A,Company_B,Company_C)
+
+# One firm
+of <- amm(amm.type="onefirm",rj=Company_A,
+            nlags=NA,
+            verbose=TRUE,
+            dates= as.Date(c("2005-01-15","2006-01-07","2007-01-06",
+                       "2008-01-05","2009-01-03")),
+           rM1=NIFTY_INDEX, others=INRUSD,
+           switch.to.innov=TRUE, rM1purge=TRUE, nlags=1)
+
+# Many firm
+mf<- amm(amm.type="manyfirm",regressand,
+            lags=NA, nlags=NA,
+            verbose=TRUE,
+            dates= as.Date(c("2005-01-15","2006-01-07","2007-01-06",
+                       "2008-01-05","2009-01-03")),periodnames=c("P1","P2","P3","P4"),
+           rM1=NIFTY_INDEX, others=INRUSD,
+           switch.to.innov=TRUE, rM1purge=TRUE)
+}
+
+\keyword{AMM}
\ No newline at end of file

Added: pkg/man/ammData.Rd
===================================================================
--- pkg/man/ammData.Rd	                        (rev 0)
+++ pkg/man/ammData.Rd	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,17 @@
+\name{ammData}
+\alias{ammData}
+
+
+\title{Three years and three random companies}
+
+\description{Data for testing and examples in the package}
+
+\usage{data(ammData)}
+
+\examples{
+library(zoo)
+data(ammData)
+str(y3c3)
+}
+
+\keyword{ammData}
\ No newline at end of file

Added: pkg/man/firmExposures.Rd
===================================================================
--- pkg/man/firmExposures.Rd	                        (rev 0)
+++ pkg/man/firmExposures.Rd	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,72 @@
+\name{firmExposures}
+\alias{firmExposures}
+
+\title{A function that estimates exposure for a single firm in a
+single period}
+
+\description{This function typically utilises an rj vector
+  and an X matrix of explanatory variables obtained 
+  using the makeX function. This would compute the exposure
+  for all columns in X.
+}
+
+\usage{
+firmExposures(rj, X, nlags = NA, verbose = FALSE)
+}
+
+\arguments{
+  \item{rj}{A vector of data for one firm}
+
+  \item{X}{A matrix of explanatory variables obtained from
+    the makeX function. The first variable is always the
+    stock market index. Other variables could be
+    such as currency or bond returns, or other variables as
+  desired by the user. 
+  }
+
+  \item{nlags}{Number of lags of explanatory variables. When unspecified
+    the best lag using the AIC is used.
+  }
+
+  \item{verbose}{Default is FALSE. When set to TRUE, the function
+    prints detailed results while running this function. 
+  }
+}
+
+
+\value{ The function returns the exposures, HAC adjusted standard
+  errors, the number of lags used, and the residuals from the fitted
+  model.
+
+\item{exposures}{This contains the exposure estimates for the firm
+j}
+
+\item{s.exposures}{This contains the HAC adjusted standard error
+of the exposures estimated for the firm rj.}
+
+\item{nlags}{Specifies the lag length provided by the user}
+
+\item{lm.res}{The model estimates used in the analysis are stored here.}
+
+}
+
+\author{Ajay Shah, Vimal Balasubramaniam}
+
+\seealso{
+\code{\link{manyfirmsAMM}}
+}
+
+\examples{
+# firmExposures
+data("firmExposuresData",package="amm")
+rj  <- test$Company_A
+rM1 <- test$NIFTY_INDEX
+rM2 <- test$usdinr
+rM3 <- test$baa
+X <- makeX(rM1, others=rM2,
+           switch.to.innov=FALSE, rM1purge=FALSE, verbose=FALSE)
+a <- firmExposures(rj, X, nlags=0, verbose=FALSE)
+print(a)
+}
+
+\keyword{firmExposures}
\ No newline at end of file

Added: pkg/man/firmExposuresData.Rd
===================================================================
--- pkg/man/firmExposuresData.Rd	                        (rev 0)
+++ pkg/man/firmExposuresData.Rd	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,19 @@
+\name{firmExposuresData}
+\alias{firmExposuresData}
+
+
+\title{Data for testing the amm package and examples therein}
+
+\description{This is a zoo data that is used for testing purposes}
+
+\usage{data(firmExposuresData)}
+
+\author{Vimal Balasubramaniam}
+
+\examples{
+library(zoo)
+data(firmExposuresData)
+str(test)
+}
+
+\keyword{firmExposuresData}
\ No newline at end of file

Added: pkg/man/makeX.Rd
===================================================================
--- pkg/man/makeX.Rd	                        (rev 0)
+++ pkg/man/makeX.Rd	2013-07-13 10:30:57 UTC (rev 73)
@@ -0,0 +1,78 @@
+\name{makeX}
+\alias{makeX}
+\title{A function to prepare explanatory variables for computation
+of Augmented Market Models
+}
+
+\description{
+  This function creates a matrix of explanatory variables in the
+  form specified by the user for further estimation of Augmented
+  market models
+}
+
+\usage{
+makeX(rM1, others, switch.to.innov = rep(TRUE, NCOL(others)), rM1purge = TRUE, nlags = 5, dates = NULL, verbose = TRUE)
+}
+
+\arguments{
+  \item{rM1}{This is generally fixed to the stock market index. The
+  first column vector of this matrix contains this variable. 
+}
+
+  \item{others}{A zoo matrix with other regressors of interest in the
+  AMM - this could be currency, bond returns, foreign flows or any
+  other variable.  
+}
+
+  \item{switch.to.innov}{is a vector of booleans with one element for
+    each column in 'others'. This controls whether or not that column in
+    'others' is switched from raw values to AR residuals. 
+}
+
+  \item{rM1purge}{whether or not the effects of all these 'others' are
+  purged from 'rM1' (i.e., is it replaced by residuals of a model
+  explaining rM1 using all these 'others'). 
+}
+
+  \item{nlags}{The number of lag terms present in this model explaining
+  rM1 using all these 'others'.
+}
+
+  \item{dates}{Specified break dates (either from structural breaks in
+  exchange rate regimes) so that all these steps are constructed within
+  each sub-period divided by the dates.
+}
+
+  \item{verbose}{Whether detailed print while running this function
+  is required.
+}
+
+}
+
+\section{Warning}{The input data has to be all na.omit(). There
+  should be no NAs in the dataset fed into this function. 
+}
+
+\value{This function generates a user specified matrix of explanatory
+  variables that will be further used in running Augmented market
+  models. 
+}
+
+
+\author{Ajay Shah, Vimal Balasubramaniam}
+
+\seealso{ \code{\link{do.one.piece}}}
+
+\examples{
+# makeX
+data("test",package="amm")
+rj  <- test$Company_A
+rM1 <- test$NIFTY_INDEX
+rM2 <- test$usdinr
+rM3 <- test$baa
+X <- makeX(rM1, others=rM2,
+           switch.to.innov=FALSE, rM1purge=FALSE, verbose=FALSE)
+}
+
+\keyword{makeX}
+\keyword{do.one.piece}
\ No newline at end of file

[TRUNCATED]

To get the complete diff run:
    svnlook diff /svnroot/eventstudies -r 73


More information about the Eventstudies-commits mailing list