install.packages("rrcov")
library(rrcov)


D <- read.table("wood.txt", header=TRUE, sep="\t")
regular_covm <- cov(D)
regular_covm
MCD_covm <- CovMcd(D,alpha=1/2) 
MCD_covm@cov


s.maha <- mahalanobis(D, center=colMeans(D), cov=regular_covm)
r.maha <- mahalanobis(D, center=MCD_covm@center, cov=MCD_covm@cov)

plot(c(1,nrow(D)),range(sqrt(c(s.maha,r.maha))),type="n",
     xlab="Observation",ylab="Mahal. distance")

points(1:nrow(D),sqrt(s.maha),col="green",pch=16)

points(1:nrow(D),sqrt(r.maha),col="blue",pch=16)


legend("topleft",col=c("green","blue"),cex=0.8,legend=c("Regular","MCD"),
       pch=c(16,16))


###Problem 2
library(mvtnorm)
set.seed(123)
n <- 200

D1 <- rmvnorm(n,mean=c(0,0),sigma=diag(2))
plot(D1, pch = 19, col="blue")

cov1  <- cov(D1)
rcov1 <- CovMcd(D1,alpha=1/2)@cov
cov1
rcov1

###Bivariate t-distribution

D2 <- rmvt(n,df=5,sigma=diag(2))
plot(D2, pch = 19, col="blue")
cov2  <- cov(D2)
rcov2 <- CovMcd(D2,alpha=1/2)@cov
cov2
rcov2

#####Bivariate Weibull-Gamma distribution

x1 <- rweibull(n,shape=1,scale=2)
x2 <- rgamma(n,shape=2,scale=1)
D3 <- cbind(x1,x2)
plot(D3, pch = 19, col="blue")

cov3 <- cov(D3)
rcov3 <- CovMcd(D3,alpha=1/2)@cov
cov3
rcov3