################################
### Function for MMA         ###
### (Mallows Model Averaging)###
################################

hansen <- function(dataframe, ycol=1, show.weights=TRUE){

# use quadprog-package
library(quadprog)

# Preparation of datamatrix
if(ycol!=1){dataframe <- cbind(dataframe[,ycol], dataframe[,-ycol]) }
avdata <- as.data.frame(dataframe)
n <- nrow(avdata)
m <- ncol(avdata)
colnames(avdata)<- paste("x",0:(m-1), sep="")
colnames(avdata)[1] <- c("y")
attach(avdata, warn.conflicts=FALSE)

### Hansens approach  ###

# let us first calculate all necessary models and the included coefficents
coeffmat <- matrix(0,ncol=m,nrow=m)
coeffmat[1,1]<-lm(avdata$y~1)$coefficients
resmat <- matrix(NA,ncol=m,nrow=n)
resmat[,1]<- lm(avdata$y~1)$residuals

for(i in 2:m){
        mycovariables <- paste("x",1:(i-1), sep="")
        myindependent <- paste("avdata$y ~")
        myformula     <- as.formula(paste(myindependent, paste(mycovariables, collapse="+") ))
        mymodel       <- lm(myformula)
        coeffmat[i,c(1:i)]  <- mymodel$coefficients
        resmat[,i]          <- mymodel$residuals
}

# now the final calculations

E <- t(resmat)%*%resmat   # this corresponds to Hansens ebar^T ebar
K <- c(seq(2,m+1))        # this corresponds to Hansens K-Vektor
mycov2 <- paste("x",1:(m-1), sep="")
myind2 <- paste("avdata$y ~")
myform2  <- as.formula(paste(myindependent, paste(mycovariables, collapse="+") ))
mymod2 <- lm(myform2)
shat <- (1/(n-(m+1)))*sum((mymod2$residuals)^2) # this corresponds to Hansens estimated Variance
Amat <- matrix(c(rep(1,m)),nrow=1,ncol=m)
Amat <- t(rbind(Amat,diag(1,m)))
bvec <-c(1,rep(0,m))
dvec <- -shat*K

weight <- solve.QP(E, dvec, Amat, bvec, meq=1)$solution
weight <- round(weight,digits=4)
weight <- matrix(c(weight),nrow=1,ncol=m)
colnames(weight)<- paste("w",0:(m-1), sep="")
rownames(weight)<- paste("Hansens weights")

estimates <- weight%*%coeffmat
colnames(estimates)<- paste("beta",0:(m-1), sep="")
rownames(estimates)<- paste("averaged est.")



# This is what the function will return
detach(avdata)
if(show.weights==TRUE){return(list(weight,estimates))}
if(show.weights==FALSE){return(estimates)}



}