LDA by hand gives different results than MASS I'm performing LDA by hand in R by following this formula:
$\delta_k(x) = x^T\Sigma^-1\mu_k-0.5*\mu_{k}^T\Sigma^-1\mu_k+log(\pi_k)$
where $\pi_k$ is the proportion of the data that's in the group $k$, $\Sigma$ is the covariance matrix that's assume to be constant through the different groups and $\mu_k$ is a vector with means of each predictor for class $k$.  
I wrote the following code to implement this formula which classify each variable based on the maximum $\delta_k(x)$ obtain
library(ISLR)
data("Default")

Y = Default$default
Xmat = Default[,-c(1:2)]
Xmat = t((t(Xmat)-colMeans(Xmat))/apply(Xmat,2,sd))

Xmat = as.matrix(Xmat)
## 1. Se estima la media de cada variable por cada variable predictora
Medias_hat = apply(Xmat,2, function(x) aggregate(x,by = list(Y),FUN = mean)[,2])
rownames(Medias_hat) = levels(Y)
## 2. Se estima la matriz de covarianza
Sigma_hat = cov(Xmat)
## 3. Se estima pi
pi_hat = table(Y)/length(Y)
## 4. Prediccion
temp_class = sapply(1:ncol(Medias_hat), function(x){
  Xmat%*%solve(Sigma_hat)%*%Medias_hat[x,]-
    as.numeric(1/2*t(Medias_hat[x,])%*%solve(Sigma_hat)%*%Medias_hat[x,])+log(pi_hat[x])
})
classificacion = apply(temp_class,1,which.max)
classificacion = ifelse(classificacion ==1,"No","Yes")

The problem i'm having is that using the package MASS I get different results and I don't understand why. The code with this package is as follow:
LdaMultiPred = lda(Y~as.matrix(Xmat))
predict(LdaMultiPred,as.data.frame(Xmat))$class

 A: So I figure out that the covariance matrix to use is not the one estimated over the hole data set but a covariance matrix that it is weighted for each class. The formula, as in The Elements of Statistical Learning is:
$\widehat{\Sigma} = \sum_{k=1}^K\sum_{g_i=k}(x_i-\widehat{\mu_k})(x_i-\widehat{\mu_k})^T/(N-K)$
I used the following R code to estimate the matrix:
temp_cov = NULL
  for(i in 1:(1+ncol(Xmat))){
    temp = NULL
    for(j in 1:length(levels(Y))){
      pos_nivel = Y == levels(Y)[j]
      if(i>ncol(Xmat)){
        temp = c(temp,(t(Xmat[pos_nivel,1]-Medias_hat[j,1])%*%(Xmat[pos_nivel,2]-Medias_hat[j,2]))
                 /(nrow(Xmat)-length(levels(Y))))
      } else {
        temp = c(temp,(t(Xmat[pos_nivel,i]-Medias_hat[j,i])%*%(Xmat[pos_nivel,i]-Medias_hat[j,i]))
                 /(nrow(Xmat)-length(levels(Y))))
      }
    }
    temp_cov = c(temp_cov,sum(temp))
  }  
  Sigma_hat = diag(temp_cov[1:ncol(Xmat)])
  Sigma_hat[lower.tri(Sigma_hat)] = temp_cov[-c(1:ncol(Xmat))]
  Sigma_hat[upper.tri(Sigma_hat)] = temp_cov[-c(1:ncol(Xmat))]

With this change the code from above yields the same result as the lda function from MASS package
