Tsukuba.Rの資料用

gistfile1.R

library(scatterplot3d)

steepest_descent <- function(epsilon,tau,x_bar1=-10,x_bar2=-20,beta=1,n=2000){
  #関数内関数の定義
  f <- function(x){x1 <- x[1];x2 <- x[2];return(100*(x1-x2)^2+(x1-1)^2)}

  #数式微分
  trig.exp <- expression(100*(x1-x2)^2+(x1-1)^2)
  D.sc1 <- D(trig.exp, "x1")
  D.sc2 <- D(trig.exp, "x2")

  #最急降下法の計算
  x_vec <- matrix(0,n+1,2,byrow=TRUE)
  for(i in seq(n)){
    nabla_f <- c((function(x1,x2){eval(D.sc1)})(x_bar1,x_bar2),
                 (function(x1,x2){eval(D.sc2)})(x_bar1,x_bar2))
    while(f(c(x_bar1,x_bar2) - beta * nabla_f) >
          f(c(x_bar1,x_bar2)) - epsilon * beta * sum(c((function(x1,x2){eval(D.sc1)})(x_bar1,x_bar2),(function(x1,x2){eval(D.sc2)})(x_bar1,x_bar2))^2)){
      beta <- tau * beta
    }
    a <- beta
    x_bar1_old <- x_bar1
    x_bar1 <- x_bar1 - a * (function(x1,x2){eval(D.sc1)})(x_bar1,x_bar2)
    x_bar2 <- x_bar2 - a * (function(x1,x2){eval(D.sc2)})(x_bar1_old,x_bar2)
    x_vec[i,] <- c(x_bar1,x_bar2)
   }

  #3次元のplot関係
  x <- seq(-1,1,length=50)
  y <- x
  z <- outer(x,y,function(x,y){mapply(function(x,y){f(c(x,y))},x,y)})


#  title(main="最急降下法",sub=paste(deparse(substitute(tau)),"=",tau,sep=""))
  cat(deparse(substitute(x)))

  persp(x,y,z,theta=150,phi=10,expand=0.5) -> res
  x <- x_vec[,1]
  y <- x_vec[,2]
  z <- mapply(function(x,y){f(c(x,y))},x,y)
  title(main="最急降下法",sub=paste(
                            deparse(substitute(tau)),"=",tau,
                            " ",
                            deparse(substitute(epsilon)),"=",epsilon,
                            sep=""))
  lines(trans3d(x,y,z, pmat = res),col="red",lwd=5)

  #結果を明示的な形で返さないようにする
  return(invisible(x_vec))
}

plot
x_vec
plot(1:10)


tau <- 0.1
substitute(tau)

deparse(substitute(tau))

png(paste(substitute(tau),"hoge.png",sep=""))

dmp <- function(x) cat(deparse(substitute(x)),"=",x,"\n")
x <- 3
dmp(x)
x = 3 
y <- "abc"
dmp(y)



dev.off()


deparse(substitute(c))
setwd("~/Desktop")

par(ask=FALSE)
mat <- matrix(1:9,3,3,byrow=TRUE)
mat
layout(mat)
for(epsilon in c(0.1,0.5,0.9)){
  for(tau in c(0.1,0.5,0.9)){
    steepest_descent(epsilon,tau)
  }
}

getwd()
















#直積を書き出す関数
write.direct.product <- function(x,y){
  e <- expand.grid(x,y)
  m <- mapply(function(x,y){paste("(",x,",",y,")",sep="")},e[,1],e[,2])
  return(paste("{",paste(m,collapse=","),"}",sep=""))
}
write.direct.product(1:3,4:6)

#R \times RからRへの写像の例
apply(expand.grid(1:3,4:6),1,function(x){x[1]*x[2]})

apply(expand.grid(letters[1:3],1:3),1,function(d){Reduce(function(x,init){paste(init,x,sep="")},rep(d[1],d[2]),init="")})

#和集合を作る操作
my.union <- function(list){
  Reduce(function(x,init){unlist(union(x,init))},list,list())
}
my.union(list(1:30,2:50,4:1,list(1,2,list(3,4))))
my.union(list(list(1,2),list(3,4)))

#共通集合を作る操作
my.intersect <- function(list){
  Reduce(function(x,init){unlist(intersect(x,init))},list,list[[1]])
}
my.intersect(list(1:30,2:50,4:10))

#冪集合を返す関数。空集合はNAにしてみた
powerset <- function(list){
  c(NA,unlist(lapply(seq(length(list)),function(x){apply(combn(list,x),2,function(y){y})}),recursive=FALSE))
}
#冪集合の例
powerset(list("a","b","c"))
lapply(powerset(list("a","b","c")),unlist)

powerset(list("a","b","c","d",1:10))

lapply(powerset(list("a","b","c","d")),unlist)[[x]]
my.union(mapply(function(x){lapply(powerset(list("a","b","c","d")),unlist)[[x]]},c(16,15)))

my.intersect(mapply(function(x){lapply(powerset(list("a","b","c","d")),unlist)[[x]]},c(16,1)))

write.list <- function(list){
  write.vector <- function(vector){
    return(paste("{",paste(vector,collapse=","),"}",sep=""))
  }
  write.list.intern <- function(list){
    ifelse(is.list(list),Reduce(function(init,x){paste(init,write.list.intern(x),sep="")},list,""),write.vector(list))
  }
  gsub("NA","",paste("{",gsub("\\}\\{","\\},\\{",write.list.intern(list)),"}",sep=""))
}

write.list(list("a","b","c","d",list(1:10,list(1,2))))
write.list(lapply(powerset(list("a","b","c","d")),unlist))


a <- "Misho"
b <-"ｷﾗｯ"
paste(a,b,sep="☆")