Machine Learning W4400 HW4

Shiyu Dong

setwd("~/Desktop/Homework_backup/machinelearning/HW4")
H <- matrix(readBin("histograms.bin", "double", 640000), 40000, 16)

H[which(H == 0, arr.in = T)] <- 0.01

MEM <- function(H, k, tau) {
    d <- ncol(H)
    n <- nrow(H)
    initial <- sample(seq(1, 40000, by = 1), k)
    t.initial <- t(apply(H[initial, ], 1, function(x) x))
    t <- matrix(nrow = k, ncol = d)
    phi <- matrix(nrow = n, ncol = k)
    a1 <- matrix(nrow = n, ncol = k)
    a <- matrix(nrow = n, ncol = k)
    c <- sample(c(1:n), k)/n
    b <- matrix(nrow = k, ncol = d)
    m <- c()
    t <- t.initial
    phi <- exp(H %*% log(t(t)))
    # for(i in 1:n){ for(kk in 1:k) { a1[i,kk]<-(c[kk]*phi[i,kk])/(phi[i,]%*%c)
    # } }
    a1 <- t(t(phi) * c) * matrix(rep(1/rowSums(phi), k), n, k)

    c <- apply(a1, 2, sum)/n
    b <- t(a1) %*% H
    t <- t(apply(b, 1, function(x) x/sum(x)))
    i = 1
    T <- c()
    for (I in 1:100) {
        N = I + 1
        phi <- exp(H %*% log(t(t)))
        # for(i in 1:n){ for(kk in 1:k) { a[i,kk]<-(c[kk]*phi[i,kk])/phi[i,]%*%c } }
        a <- t(t(phi) * c) * matrix(rep(1/rowSums(phi), k), n, k)

        c <- apply(a, 2, sum)/n
        b <- t(a) %*% H
        t <- t(apply(b, 1, function(x) x/sum(x)))
        Tau <- norm(a - a1)
        T[I] <- Tau
        if (Tau < tau) {
            m <- apply(a, 1, which.max)
            result <- list(M = m, N, TAU = T)
            return(result)
        }
        a1 <- a
    }
}

k=3

par(mfrow = c(2, 2))
tau <- c(0.001, 1e-04, 1e-05, 1e-06)
k <- 3
for (t in tau) {
    m <- MEM(H, k, t)
    ima <- matrix(m[[1]], nrow = 200, byrow = T)
    mm <- cbind(seq(1, 200, by = 1), ima)
    image <- mm[order(-mm[, 1]), ][, -1]
    title = paste("k=3 tau=", t)
    image(x = 1:200, y = 1:200, t(image), axes = FALSE, col = grey((0:256)/256), 
        xlab = "", ylab = "", main = title)
}

k=4

par(mfrow = c(2, 2))
tau <- c(0.001, 1e-04, 1e-05, 1e-06)
k <- 4
for (t in tau) {
    m <- MEM(H, k, t)
    ima <- matrix(m[[1]], nrow = 200, byrow = T)
    mm <- cbind(seq(1, 200, by = 1), ima)
    image <- mm[order(-mm[, 1]), ][, -1]
    title = paste("k=4 tau=", t)
    image(x = 1:200, y = 1:200, t(image), axes = FALSE, col = grey((0:256)/256), 
        xlab = "", ylab = "", main = title)
}

k=5

par(mfrow = c(2, 2))
tau <- c(0.001, 1e-04, 1e-05, 1e-06)
k <- 5
for (t in tau) {
    m <- MEM(H, k, t)
    ima <- matrix(m[[1]], nrow = 200, byrow = T)
    mm <- cbind(seq(1, 200, by = 1), ima)
    image <- mm[order(-mm[, 1]), ][, -1]
    title = paste("k=5 tau=", t)
    image(x = 1:200, y = 1:200, t(image), axes = FALSE, col = grey((0:256)/256), 
        xlab = "", ylab = "", main = title)
}