#################################################
# A04: R程式設計 (& R Style Guide)              #
# 吳漢銘 國立政治大學統計學系                   #
# https://hmwu.idv.tw                           #
#################################################

# 3/110
getwd()
dir()
# setwd("C:\\Program Files\\R\\working")


# 5/110
x <- 1
if((x-2) < 0) cat("expr1 \n") else cat("expr2 \n")
if((x-2) > 0) cat("expr1 \n") else cat("expr2 \n")

x <- c(-1, 2, 3)
if((x-2) < 0) cat("expr1 \n") else cat("expr2 \n")
if((x-2) > 0) cat("expr1 \n") else cat("expr2 \n")


# 6/110
x <- 0
if(x) cat("expr1 \n") else cat("expr2 \n")
if(x+1) cat("expr1 \n") else cat("expr2 \n")

x <- c(-1, 0, 1, 2, 3)
if(x) cat("expr1 \n") else cat("expr2 \n")
if(x+1) cat("expr1 \n") else cat("expr2 \n")


# 7/110
x <- c(-1, 2, 3)
if(any(x <= 0)) y <- log(1+x) else y <- log(x)
y
z <- if(any(x <= 0)) log(1+x) else log(x) 
z

x <- c(-1, 2, 3)
if(any(x <= 0)){ 
     y <- log(1+x) 
 } else{ 
    y <- log(x)
 }
y

x <- c(-1, 2, 3)
if(any(x <= 0)){ 
  y <- log(1+x) 
} 
else{ # wrong code for demo
  y <- log(x)
}
y


# 8/110
(a <- sample(1:5))
(b <- sample(1:5))
a == b
any(a == b)
all(a == b)


check.if <- function(a, b){
  if(a == b){
    cat("Equal! \n")
  }else{
    cat("Not equal! \n")  
  }
}

check.if2 <- function(a, b){
  if(sum(abs(a - b)) == 0){
    cat("Equal! \n")
  }else{
    cat("Not equal! \n")  
  }
}

a <- sample(1:10, 4)
b <- a
check.if2(a, b)

check.if(a = 1, b = 1)
check.if(a = 1, b = 2)
check.if(a = 1, b = c(1, 2, 3))
check.if(a = 1, b = c(2, 1, 3))
check.if(a = c(1, 2, 3), b = c(1, 2, 3))
check.if(a = c(2, 4, 5), b = c(1, 2, 3))
check.if(a = c(1, 5), b = c(4, 2, 3))

identical(a, b) 
all.equal(pi, 355/113) 


# 10/110
x <- 3
y <- 4
x < 2
y > 2
x < 2 || y > 2
x > 2
y > x
x > 2 && y > x

x < 2 | y > 2
x > 2 & y > x

xv <- c(1, 2, 3)
yv <- c(2, 2, 5)
xv < 2
yv > 2
xv < 2 || yv > 2
(! xv < 2) || yv > 2
xv < 2 || (! yv > 2)
xv < 2 && yv > 2
(! xv < 2) && yv > 2
xv < 2 && (! yv > 2)

xv < 2 | yv > 2
(! xv < 2) | yv > 2
xv < 2 | (! yv > 2)
xv < 2 & yv > 2
(! xv < 2) & yv > 2
xv < 2 & (! yv > 2)


# 11/110
a <- 2.13

if( a > 10 ){ 
    cat("a > 10 \n")
}else if(a > 5){ 
    cat("5 < a < 10 \n")
}else if(a > 2.5){ 
    cat("2.5 < a < 5 \n")
}else if(a > 1.25){
    cat("1.25 < a < 2.5 \n")
}else{
    cat("a < 1.25")
}

1.25 < a < 2.5 

x <- 0
if (x < 0) {
   print("Negative number")
} else if (x > 0) {
   print("Positive number")
} else
   print("Zero")


# 12/110
(x <- c(2:-1))
sqrt(x)
sqrt(ifelse(x >= 0, x, NA))
ifelse(x >=  0, sqrt(x), NA)

(yes <- 5:6)
(no <- pi^(0:2))
ifelse(NA, yes, no)
ifelse(TRUE, yes, no)
ifelse(FALSE, yes, no)


# 13/110
x <- c(24, 13, 26, 21,  7,  9, 2, 1, 30, 14, 20, 16,  
       6, 4, 12, 8, 11, 22, 18, 3)
ifelse(x <= 10, 1, ifelse(x <= 20, 2, 3))

set.seed(12345)
age <- sample(1:100, 20)
age

set.seed(12345)
code <- sample(LETTERS[1:5], 20, replace = T)
code


# 14/110
set.seed(12345)
score <- sample(0:100, 10, replace = T)
score

gpa.table <- data.frame(grade = c("A", "B", "C", "D", "E"), 
                        pscore = c("80-100", "70-79", "60-69", "50-59", "49-0"), 
                        GPA = c(4, 3, 2, 1, 0))
gpa.table
set.seed(12345)
score <- sample(0:100, 10, replace = T)

score_to_gpa <- function(x){

    group.id <- ifelse(x >=  80, 1, 
                   ifelse(x >=  70, 2, 
                          ifelse(x >= 60, 3, 
                                 ifelse(x >=  50, 4, 5))))
    data.frame(score = x, gpa.table[group.id, ], row.names = NULL)
}

score_to_gpa(score)


# 20/110
my_dist <- function(x1, y1, x2, y2){
  d <- sqrt((x1-x2)^2 + (y1-y2)^2)
  d
}

my_dist2 <- function(x1, y1, x2 = 0, y2 = 0){
  d <- sqrt((x1-x2)^2 + (y1-y2)^2)
  list(points.a = c(x1, y1), points.b = c(x2, y2), dist.ab = d)
}

my_dist(1, 2, 4, 7)

my_dist2(1, 2, 4, 7) 

my_dist2(1, 2) 


# 21/110
my_dist <- function(a, b){
   sqrt(sum((a-b)^2))
}

my_dist(a = c(1, 2), b = c(4, 7))

f <- function(x){
  x^2+1
}
 
x <- 1:5
y <- f(x)
y

x <- 1:5
y <- x^2+1
y


# 22/110
min(5:1, pi)
pmin(5:1, pi)

parmax <- function(a, b){
	c <- pmax(a, b)
	median(c)
}
x <- c(1, 9, 2, 8, 3, 7)
y <- c(9, 2, 8, 3, 7, 2)
parmax(x, y)

data_ratio <- function(x){
    x.number <- length(x)
    x.up <- mean(x) + sd(x)
    x.down <- mean(x) - sd(x)
    x.n <- length(x[x.down < x & x < x.up])
    x.p <- x.n/x.number
    list(number = x.n, percent = x.p)
}


data_ratio(iris[, 1])


# 23/110
compute <- function(a, b = 0.5){

    sum <- a + b
    diff <- a - b
    prod <- a * b
    if(b != 0){
      div <- a / b
    }else{
      div <- "divided by zero"
    }
    list(sum = sum, diff = diff, product = prod, divide = div)    
}

norm <- function(x) sqrt(x%*%x)
norm(1:4)

compute(2, 5)

compute(2)

compute(2, 0)


# 24/110
two_sample_test <- function(y1, y2){
    n1 <- length(y1); n2 <- length(y2)
    m1 <- mean(y1); m2 <- mean(y2)
    s1 <- var(y1); s2 <- var(y2)
    s <- ((n1-1)*s1 + (n2-1)*s2)/(n1+n2-2)
    stat <- (m1-m2)/sqrt(s*(1/n1+1/n2))
    list(means = c(m1, m2), pool.var = s, stat = stat)
}
t.stat <- two_sample_test(iris[, 1], iris[, 2]) 
t.stat 
t.stat 


# 25/110
rm(list = ls())
my_sqrt_sum <- function(x, y){
    a <- sqrt(x)
    b <- sqrt(y)
    c <- a+b
    c
}

a <- 4
b <- 9
my_sqrt_sum(a, b)
a
b

rm(list = ls())
my_sqrt_sum <- function(x, y){
    a <- sqrt(x)
    b <- sqrt(y)
    c <- a+b
    c
}
my_sqrt_sum(4, 9)
a
b


# 26/110
rm(list = ls())
y <- 9
my_sqrt_sum <- function(x){
    a <- sqrt(x) 
    b <- sqrt(y)
    y <- sqrt(y)
    c <- a+b
    c
}
my_sqrt_sum(4)
a
b
y

rm(list = ls())
Y.VALUE <- 9
my_sqrt_sum <- function(x){
    a <- sqrt(x) 
    b <- sqrt(Y.VALUE)
    c <- a+b
    c
}
my_sqrt_sum(4)

rm(list = ls())
my_sqrt_sum <- function(x, y){
    x <- sqrt(x)
    y <- sqrt(y)
    c <- x+y
    c
}

x <- 4
y <- 9
x <- my_sqrt_sum(x, y)
x
y


# 27/110
myfun1 <- function(x){   
   y <- x + 5
   cat("y: ", y, "\n")
}

myfun2 <- function(x){   
   y <<- x + 5
   cat("y: ", y, "\n")
}

y <- 5; cat("y: ", y, "\n")
myfun1(3)
cat("y: ", y, "\n")

y <- 5; cat("y: ", y, "\n")
myfun2(3)
cat("y: ", y, "\n")

myfun1 <- function(x){   
   y <- x + 5
   cat("y: ", y, "\n")
}
 
myfun2 <- function(x){   
   y <<- x + 5
   cat("y: ", y, "\n")
}
 
y <- 5; cat("y: ", y, "\n")
myfun1(3)
cat("y: ", y, "\n")
 

y <- 5; cat("y: ", y, "\n")
myfun2(3)
cat("y: ", y, "\n")


# 28/110
my_stat <- function(x){
    x.number <- length(x)
    x.mean <- mean(x)
    x.sd <- sd(x)
    list(number = x.number, mean = x.mean, sd = x.sd)
}

my_stat(iris[, 1])


# 29/110
lm

myfun <- function(x, ...){   
    y <- mean(...) + x
    y
}
 
data <- rnorm(40)
myfun(6, data)


# 30/110
many_means <- function(...){

    #use [[]] subscripts in addressing its elements.
    data <- list(...) 
    n <- length(data)
    means <- numeric(n)
    vars  <- numeric(n)
    for(i in 1:n){
        means[i] <-  mean(data[[i]])
	vars[i] <-  var(data[[i]])
    }

    print(means)
    print(vars)
}

x <- rnorm(100); y <- rnorm(200); z <- rnorm(300)
many_means(x, y, z)


# 31/110
data_k_ratio <- function(x, k = 1){
   x.number <- length(x)
   x.mean <- mean(x)
   x.sd <- sd(x)
   x.up <- x.mean + k*x.sd;
   x.down <- x.mean - k*x.sd;
   x.n <- length(x[(x.down < x) & (x < x.up)])
   x.p <- x.n/x.number
   list(number = x.n, percent = x.p)
}
library(MASS)
data_k_ratio(drivers, 1)
data_k_ratio(drivers, 2)
data_k_ratio(drivers, 3)

library(MASS)
data_k_ratio(drivers, 1)

data_k_ratio(drivers, 2)

data_k_ratio(drivers, 3)


# 33/110
for(i in 1:5){
    cat("loop: ", i, "\n")
}

for(k in c(1, 17, 3, 56, 2)){
    cat(k, "\t")
}

for(bloodType in c("A", "AB", "B", "O")){
    cat(bloodType, "\t")
}


# 34/110
rm(list = ls())
y <- round(rnorm(10), 2)
z <- y
y
i
for(i in 1:length(y)){
	if(y[i] < 0) 
	y[i] <-0
}
y
i

z[z<0] <- 0
z

rm(list = ls())
y <- round(rnorm(10), 2)
z <- y
y
i
for(i in 1:length(y)){
if(y[i] < 0) 
y[i] <-0
}
y
i

z[z < 0] <- 0
z


# 35/110
a <- numeric(5)
for(i in 1:5){
  a[i]<- i^2
}
 a


m <- 3
n <- 4
for(i in 1:m){
   for(j in 1:n){
     cat("loop: (", i, ", ", j, ")\n")
   }
}

a <- matrix(0, 2, 4)
for(i in 1:2){
    for(j in 1:4){ 
       a[i, j]<- i+j
    }
}
a


# 36/110
m <- 3
n <- 4
for(i in 1:m){
   for(j in 1:n){

     if(i == 2){
       cat("before next:", i, ", ", j, "\n")
       next
       cat("after next:", i, ", ", j, "\n")
     }else{
       cat("loop: (", i, ", ", j, ")\n")
     }
   }
}


# 37/110
m <- 3
n <- 4
for(i in 1:m){
   for(j in 1:n){

     if(i == 2){
       cat("before break:", i, ", ", j, "\n")
       break
       cat("after break:", i, ", ", j, "\n")
     }else{
       cat("loop: (", i, ", ", j, ")\n")
     }
   }
}


# 38/110
check_prime <- function(num){
  
  yes <- FALSE
  
  if(num == 2){
    yes <- TRUE

  } else if(num > 2) {
    
    yes <- TRUE
    for(i in 2:(num-1)) {
        if ((num %% i) == 0) {
            yes <- FALSE
            break
        }
    }
  } 
  
  if(yes) {
    cat(num, "is a prime number. \n")
  } else {
    cat(num, "is not a prime number. \n")
  }
}

check_prime(2)
check_prime(13)
check_prime(25)


# 40/110
a <- 5
while(a > 0){  
    a <- a-1  
    cat(a, "\n")
    if(a == 2){
        cat("before next:", a, "\n")
        next
        cat("after next:", a, "\n")
    }
}

a <- 5
while(a > 0){  

    if(a == 2){
        cat("before break:", a, "\n")
        break
    }
    a <- a-1  
    cat(a, "\n")
}

a <- 5
while(a > 0){  

    if(a == 2){
        cat("before break:", a, "\n")
        next
        cat("after break:", a, "\n")
    }
    a <- a-1  
    cat(a, "\n")
}


# 41/110
factorial_for <- function(n){
    f <- 1
    if(n < 2) return(1)
    for(i in 2:n){
        f <- f * i
    }
    f
}
factorial_for(5)

factorial_while <- function(n){
    f <- 1
    t <- n
    while(t > 1){
        f <- f * t
        t <- t - 1
    }
    return(f)
}
factorial_while(5)

factorial_repeat <- function(n){
    f <- 1
    t <- n
    repeat{
        if(t < 2) break
        f <- f * t
        t <- t - 1
    }
    return(f)
}
factorial_repeat(5)

factorial_call <- function(n, f){     
    
    if(n <=  1){
      return(f)
    }
    else{ 
      factorial_call(n - 1, n * f)
    }
}
factorial_call(5, 1)

factorial_cumprod <- function(n) max(cumprod(1:n))
factorial_cumprod(5)
factorial(5)


# 42/110
x <- 3
switch(x, 
         cat("2+2\n"), 
         cat("mean(1:10)\n"), 
         cat("sd(1:10)\n"))

switch(x, 2+2, mean(1:10), sd(1:10))

switch(2, 2+2, mean(1:10), sd(1:10))

switch(6, 2+2, mean(1:10), sd(1:10))


my_lunch <- function(y){
    switch(y, 
           fruit = "banana", 
	   vegetable = "broccoli", 
	   meat = "beef")
}
my_lunch("fruit")
my_lunch(fruit)


# 43/110
x_center <-  function(x, type){
    switch(type, 
           mean = mean(x), 
           median = median(x), 
           trimmed = mean(x, trim = 0.1), 
           stop("Measure is not included!"))
}

x <- rnorm(20)
x_center(x, "mean")
x_center(x, "median")
x_center(x, "trimmed")
x_center(x, "mode")



# 44/110
my_median_1 <- function(x){
    odd.even <- length(x)%%2
    if(odd.even == 0){
         (sort(x)[length(x)/2] + sort(x)[1+length(x)/2])/2
    }else{
         sort(x)[ceiling(length(x)/2)]
    }
}

my_median_2 <- function(x){
    odd.even <- length(x)%%2
    s.x <- sort(x)
    n <- length(x)
    if(odd.even == 0){
        median <- (s.x[n/2] + s.x[1+n/2])/2
    }else{
        median <- s.x[ceiling(n/2)]
    }    
    return(median)
}

x <- rnorm(30)
my_median_1(x) 
my_median_2(x) 
median(x)


# 45/110
(x <- matrix(1:24, nrow = 4))
apply(x, 1, sum)  
apply(x, 2, sum)
apply(x, 1, sqrt) 
apply(x, 2, sqrt)


# 46/110
math <- sample(1:100, 50, replace = T)
english <- sample(1:100, 50, replace = T)
algebra <- sample(1:100, 50, replace = T)
ScoreData <- cbind(math, english, algebra)
head(ScoreData, 5)
myfun <- function(x){
     sqrt(x)*10
}
sdata1 <- apply(ScoreData, 2, myfun)
head(sdata1, 5)

head(apply(ScoreData, 2, function(x) sqrt(x)*10), 5)
myfun2 <- function(x, attend){
    y <- sqrt(x)*10 + attend
    ifelse(y > 100, 100, y)
}
sdata2 <- apply(ScoreData, 2, myfun2, attend = 5)
head(sdata2, 5)


# 47/110
tapply(iris$Sepal.Width, iris$Species, mean)
set.seed(12345)
scores <- sample(0:100, 50, replace = T)
grade <- as.factor(sample(c("大一", "大二", "大三", "大四"), 50, replace = T))
bloodtype <- as.factor(sample(c("A", "AB", "B", "O"), 50, replace = T))
tapply(scores, grade, mean)
tapply(scores, bloodtype, mean)
tapply(scores, list(grade, bloodtype), mean)


# 48/110
n <- 20
(my.factor <- factor(rep(1:3, length = n), levels = 1:5))
table(my.factor)
tapply(1:n, my.factor, sum)

tapply(1:n, my.factor, range)

tapply(1:n, my.factor, quantile)

# not run
# > by(iris[, 1:4] , iris$Species , mean)
by(iris[, 1:4] , iris$Species , colMeans)

varMean <- function(x, ...) sapply(x, mean, ...)
by(iris[, 1:4], iris$Species, varMean)


# 49/110
a <- c("a", "b", "c", "d")
b <- c(1, 2, 3, 4, 4, 3, 2, 1)
c <- c(T, T, F)
 
list.object <- list(a, b, c)
 
my.la1 <- lapply(list.object, length)
my.la1

my.la2 <- lapply(list.object, class)
my.la2


# 50/110
rep(5.6, 3)
replicate(3, 5.6)
rep(rnorm(1), 3)
replicate(3, rnorm(1))
replicate(3, mean(rnorm(10)))

dice1 <- sample(1:6, 1)
dice2 <- sample(1:6, 1)
dice1 + dice2

my_dice <- function(n){
   dice.no <- sample(1:6, n)
   dice.sum <- sum(dice.no)
   output <- c(dice.no, dice.sum)
   names(output) <- c(paste0("dice", 1:n), "sum")
   output
 }
my_dice(3)
replicate(5, my_dice(2))


# 52/110
(select.num <- sapply(iris, is.numeric))  
iris[1:2, select.num]
select.fac <- sapply(iris, is.factor)
select.fac
iris[1:5, select.fac]
apply(iris, 2, is.numeric)
unique(iris$Species)
table(iris$Species)


# 53/110
wf <- read.table("worldfloras.txt", header = TRUE)
attach(wf)
names(wf)
dim(wf)


# 54/110
index <- grep("R", as.character(Country))
as.vector(Country[index])
as.vector(Country[grep("^R", as.character(Country))])
as.vector(Country[grep(" R", as.character(Country))])
as.vector(Country[grep(" ", as.character(Country))])
as.vector(Country[grep("y$", as.character(Country))])


# 55/110
my.pattern <- "[C-E]" 
index <- grep(my.pattern, as.character(Country))
as.vector(Country[index]) 

as.vector(Country[grep("^[C-E]", as.character(Country))])

as.vector(Country[-grep("[a-t]$", as.character(Country))])

as.vector(Country[-grep("[A-T a-t]$", as.character(Country))])


# 56/110
as.vector(Country[grep("^.y", as.character(Country))])

as.vector(Country[grep("^..y", as.character(Country))])

as.vector(Country[grep("^.{5}y", as.character(Country))])

as.vector(Country[grep("^.{, 4}$", as.character(Country))]) 

as.vector(Country[grep("^.{15, }$", as.character(Country))])


# 57/110
text <- c("arm", "leg", "head", "foot", "hand", "hindleg", "elbow")
text
gsub("h", "H", text)
gsub("o", "O", text)
sub("o", "O", text)
gsub("^.", "O", text)


x <- c(3, 2, 1, 0, 4, 0)
replace(x, x == 0, 1)

replace(text, text == "leg", "LEG")
replace(text, text %in% c("leg", "foot"), "LEG")


# 58/110
text <- c("arm", "leg", "head", "foot", "hand", "hindleg", "elbow")
regexpr("o", text)
grep("o", text)
text[grep("o", text)]
gregexpr("o", text)

charmatch("m", c("mean", "median", "mode"))
charmatch("med", c("mean", "median", "mode"))


# 59/110
stock <- c("car", "van")
requests <- c("truck", "suv", "van", "sports", "car", "waggon", "car")
requests %in% stock
index <- which(requests %in% stock)
requests[index]

x <- round(rnorm(10), 2)
x
index <- which(x < 0)
index
x[index]
x[x < 0]


# 60/110
x <- matrix(sample(1:12), ncol = 4, nrow = 3)
x
which(x %% 3 == 0)
which(x %% 3 == 0, arr.ind = T)

x <- c(45, 3, 50, 41, 14, 50, 3)
which.min(x)
which.max(x)
x[which.min(x)]
x[which.max(x)]
which(x == max(x))

match(1:10, 4)
match(1:10, c(4, 2))
x
match(x, c(50, 3))


# 61/110
setA <- c("a", "b", "c", "d", "e")
setB <- c("d", "e", "f", "g")

union(setA, setB)
intersect(setA, setB)
setdiff(setA, setB)
setdiff(setB, setA)
setA %in% setB
setB %in% setA
setA[setA %in% setB] 


# 62/110
myFun <- function(n){
  x <- 0
  for(i in 1:n){
    x <- x + i
  }
  x
}

system.time({
     ans <- myFun(10000)
 })

start.time <- proc.time()
for(i in 1:50) mad(runif(500))
proc.time() - start.time

start.time <- Sys.time()
ans <- myFun(10000)
end.time <- Sys.time()
end.time - start.time


#install.packages("profvis")
library(profvis)

profvis({
  
  data(diamonds, package = "ggplot2")
  
  plot(price ~ carat, data = diamonds)
  m <- lm(price ~ carat, data = diamonds)
  abline(m, col = "red")
})


# 63/110
city <- read.table("city.txt", header = TRUE, row.names = NULL, sep = "\t")
attach(city)
names(city)

rank.price <- rank(price)
sorted.price <- sort(price)
ordered.price <- order(price)

sort(price, decreasing = TRUE)
rev(sort(price))


# 64/110
city

(view1 <- data.frame(location, price, rank.price))

(view2 <- data.frame(sorted.price, ordered.price)) 

(view3 <- data.frame(location[ordered.price], price[ordered.price])) 

# 65/110
y <- 1:20
sample(y)
sample(y)
sample(y, 5)
sample(y, 5)
sample(y, 5, replace = T)

substr("this is a test", start = 1, stop = 4)
substr(rep("abcdef", 4), 1:4, 4:5)

x <- c("asfef", "qwerty", "yuiop[", "b", "stuff.blah.yech")
substr(x, 2, 5)
substring(x, 2, 4:6)
substring(x, 2) <- c("..", "+++")
x


# 66/110
Triangular <- function(u){
	s <- ifelse(abs(u) <=  1, 1, 0)
	ans <- (1-abs(u))*s
	return(ans)
}
Gaussian <- function(u){
	ans <- exp((-1/2)*(u^2))/sqrt(2*pi)
	return(ans)
}
Epanechnikov <- function(u){
	s <- ifelse(abs(u) <=  sqrt(5) , 1, 0)
	ans <- 3*(1-((u^2)/5))/(4*sqrt(5))*s
	return(ans)
}

par(mfrow = c(1, 3))
x <- seq(-3, 3, 0.1)
plot(x, Triangular(x), main = "Triangular Kernel", type = "l")
plot(x, Gaussian(x), main = "Gaussian Kernel", type = "l")
plot(x, Epanechnikov(x), main = "Epanechnikov Kernel", type = "l")


# 67/110
fh <- function(xi, x, h, kernel, n = 150){
	ans <- sum(kernel((x-xi)/h))/(n*h)
	return(ans)
}

x <- iris[, 1]
fh(x, 7, 0.2736, Triangular)
fh(x, 7, 0.2736, Gaussian)
fh(x, 7, 0.2736, Epanechnikov)


# 68/110
binomial <- function(k, n, p){
  factorial(n)/(factorial(k) * factorial(n - k)) * (p^k) * ((1-p)^(n-k))    
}

compute_mu_sigma <- function(pmf, parameter){

  mu <- 0
  sigma2 <- 1
  
  pmf.name <- deparse(substitute(pmf))
  cat("Input is", pmf.name, "distribution.\n")
  if(pmf.name ==  "binomial"){
    k <- parameter[[1]]
    n <- parameter[[2]]
    p <- parameter[[3]]
    mu <- sum(k * pmf(k, n, p))
    sigma2 <- sum((k - mu)^2 * pmf(k, n, p))
  }  
  cat("mu: ", mu, "\n")
  cat("sigma2: ", sigma2, "\n")
}

compute_mu_sigma(pmf = binomial, parameter = c(4, 10, 0.5))


# 69/110
binomial <- function(k, n, p){
  factorial(n)/(factorial(k) * factorial(n - k)) * (p^k) * ((1-p)^(n-k))    
}

poisson <- function(k, lambda){
  exp(-lambda) * (lambda^k)/(factorial(k))   
}

geometric <- function(k, p){
  (1 - p)^k * p  
}

compute_mu_sigma <- function(pmf, parameter){
  pmf.name <- deparse(substitute(pmf))
  mu <- sum(parameter$k * (do.call("pmf", parameter)))
  sigma2 <- sum((parameter$k - mu)^2 * do.call("pmf", parameter))
  cat("distribution:　", pmf.name, "\n")
  cat("mu: ", mu, "\t sigma2:", sigma2, "\n" )
}

my.par <- list(k = c(0:10), n = 10, p = 0.6)
compute_mu_sigma(pmf = binomial, parameter = my.par)

my.par <- list(k = c(0:100), lambda = 4)
compute_mu_sigma(pmf = poisson, parameter = my.par)

my.par <- list(k = c(0:10000), p = 0.4)
compute_mu_sigma(pmf = geometric, parameter = my.par)


# 70/110
as.numeric(factor(c("a", "b", "c")))
as.numeric(c("a", "b", "c")) #don’t work


# 71/110
(x <- sample(1:42, 6))
(y <- letters)
get("x")
get("y")[1:5]
 
for(i in 1:5){
     x.name <- paste("x", i, sep = ".")
     assign(x.name, 1:i)
     cat(x.name, ": \t")
     cat(get(x.name), "\n")
 }

a <- 100
(my.math <- c("3+4", "a/5"))
eval(my.math)
eval(parse(text = my.math[1]))
plot.type <- c("plot", "hist", "boxplot")
x <- rnorm(100)
my.plot <- paste(plot.type, "(x)", sep = "")
eval(parse(text = my.plot))


# 72/110
library(e1071)
fclustIndex
plot
methods(plot)
plot.lm


# 73/110
plot.table
?plot.table
graphics:::plot.table

anova
methods(anova)
stats:::anova.nls
stats:::anova.loess

svm
?svm
methods(svm)
e1071:::svm.default

princomp 
methods(princomp) 
getAnywhere('princomp.default') # or
stats:::princomp.default


# 99/110
cat("第一題")
string <- "((1+2)*(3+4)*(5+6))/(7+8)" 
gregexpr("[(]", string)[[1]]
length(gregexpr("[(]", string)[[1]])


# 100/110
cat("第一題\n")
string <- "((1+2)*(3+4)*(5+6))/(7+8)" 
left.num <- length(gregexpr("[(]", string)[[1]])
cat("left.num: ", left.num, "\n")
right.num <- length(gregexpr("[)]", string)[[1]])
cat("right.num: ", right.num, "\n")
if(left.num ==  right.num){
   cat("OK")
} else{
   cat("Not OK")
}


# 101/110
ex1 <- function(){
    cat("第一題\n")

    #string <- "((1+2)*(3+4)*(5+6))/(7+8)" 
    cat("輸入包含左右小括號之字串(最長為40字元)，請判斷是否左右小括號配對正確")    
    string <- scan(what = "character", nmax = 1, quiet = TRUE)

    left.num <- length(gregexpr("[(]", string)[[1]])
    #cat("left.num: ", left.num, "\n")
    right.num <- length(gregexpr("[)]", string)[[1]])
    #cat("right.num: ", right.num, "\n")

    if(left.num ==  right.num){
       cat("OK")
    }else{
       cat("Not OK")
    }
}

ex1()

# 102/110
ex1 <- function(){
   cat("第一題\n")
   cat("########################################\n")
   cat("# 輸入包含左右小括號之字串(最長為40字元)，   #\n")
   cat("# 請判斷是否左右小括號配對正確              #\n")
   cat("# 例如輸入： {\tt ((1+2)-3)*(4/5)}       #\n")
   cat("########################################\n")

   ##輸入
   string <- scan(what = "character", nmax = 1, quiet = TRUE)

   ##找出"(" ")"，並計數
   left.num <- length(gregexpr("[(]", string)[[1]])
   right.num <- length(gregexpr("[)]", string)[[1]])

   ##判斷是否相等
   if(left.num ==  right.num){
     ##是的話，輸出OK
     cat("OK")
   }
   else{
     ##不是的話，輸出NOT OK
     cat("Not OK")
   }
}
Y.or.N <- "y"

while(Y.or.N == "y"){
   ex1()
   cat("繼續練習這一題(Y/N): ")
   Y.or.N <- scan(what = "character", nmax = 1, quiet = TRUE)
   if(Y.or.N != "y" & Y.or.N != "n"){
      cat("輸入錯誤，再輸入一次 ")
   }
}


# 104/110
source("http://www.hmwu.idv.tw/web/R/Rscript/ex1.R")
# source("ex1.R")