lec03

Size: px

Start display at page:

Download "lec03"

よいかずほうねん
6 years ago
Views:

2 # Previous inappropriate version tol = 1e-7; grad = 1e10; lambda = 0.2; gamma = 0.9 x = 10; x.hist = x; v = 0 while (abs(grad)>tol){ grad = 2*x+2 v = gamma*v-lambda*grad x = x + v x.hist=c(x.hist,x) print(c(x,grad)) x.temp=seq(-10,10,0.1) plot(x.temp, x.temp^2+2*x.temp,type='l',lwd=2) lines(x.hist,x.hist^2+2*x.hist,type='o',pch=1,col='red',cex=1) [1] [1] [1] e e-15

3 # New appropriate version tol = 1e-7; grad = 1e10; lambda = 0.2; gamma = 0.9 x = 10; x.hist = x; v = 0 repeat { grad = 2*x+2 if (abs(grad) <= tol){break v = gamma*v-lambda*grad x = x + v x.hist=c(x.hist,x) print(c(x,grad)) x.temp=seq(-10,10,0.1) plot(x.temp, x.temp^2+2*x.temp,type='l',lwd=2) lines(x.hist,x.hist^2+2*x.hist,type='o',pch=1,col='red',cex=1) [1] [1]

5 naive.stochoptim <- function(obj.func, x, n.iter, width){ opt.value <- do.call(obj.func, list(x)) opt.hist = matrix(0, nrow = n.iter, ncol = 5) opt.hist[1,] = c(x, x, opt.value, opt.value, 1) for (i.iter in 2:n.iter){ accpet = 0 temp.x <- x + rnorm(1, mean = 0, sd = width) temp.value <- do.call(obj.func, list(temp.x)) if (temp.value < opt.value){ x = temp.x opt.value = temp.value accept = 1 opt.hist[i.iter, ] = c(x, temp.x, opt.value, temp.value, 1) return(data.frame(opt.hist))

6 set.seed(50) n.iter =500 fun01<-function(x){x^2+2*x res <- naive.stochoptim(fun01,3,n.iter,1) > head(res) X1 X2 X3 X4 X

8 x "#$ = x + g Δx. T Δf = f x "#$ f x P x "#$, T = T ηt 0# T: temperature c: problem specific constant P: probability of accepting hypothetical x Eta: annealing constant (<1)

9 Δf = f x "#$ P x "#$, T = f x exp Δf ct f x "#$ x "#$ f x "#$ x "#$ x "#$ x "#$

10 SimAneal01<-function(func, initxy, maxitr=1000, C=1, eta=0.99, width=10) { x=initxy opt.value = do.call(func,list(x)) n.var = length(x) opt.hist=matrix(0, nrow=maxitr, ncol=5) opt.hist[1,]=c(x,x,opt.value,opt.value,0) for (i_loop in 2:maxItr) { accept = 0 temp.x = x + rnorm(n.var, mean = 0, sd=width) temp.value= do.call(func,list(temp.x)) delta=temp.value-opt.value; prob=1(1+exp(delta(c*width))) if (runif(1) < prob) { x = temp.x; opt.value = temp.value; accept = 1 opt.hist[i_loop,]=c(x, temp.x, opt.value, temp.value, accept); width=width*eta return(data.frame(opt.hist))

11 set.seed(50) n.iter =500 fun01<-function(x){x^2+2*x res <- SimAneal01(fun01, 3, n.iter, 1, 0.985, 1) > head(res) X1 X2 X3 X4 X

14 model IV1 IV2 IV3 IV4 IV5 IV6 IV7 IV8 IV9 IV

15 set.seed(19) nvar = 10 ndata = 50 # X=matrix(rnorm(nVar*nData,mean=10,sd=2),ncol=nVar); # y=x%*%c(-2,0,2,0,2,0,-3,0,-2,0)+rnorm(ndata,mean=0,sd=4); # memo # %*%

16 > summary(lm(y~x)) Coefficients: Estimate Std. Error t value Pr(> t ) (Intercept) X e-13 *** X X e-12 *** X X e-12 *** X X < 2e-16 *** X X e-15 *** X Signif. codes: 0 *** ** 0.01 * Residual standard error: on 39 degrees of freedom Multiple R-squared: , Adjusted R-squared: F-statistic: on 10 and 39 DF, p-value: < 2.2e-16

17 > summary(lm(y~x[,seq(1,9,2)])) Coefficients: Estimate Std. Error t value Pr(> t ) (Intercept) X[, seq(1, 9, 2)] e-16 *** X[, seq(1, 9, 2)] e-14 *** X[, seq(1, 9, 2)] e-14 *** X[, seq(1, 9, 2)] < 2e-16 *** X[, seq(1, 9, 2)] < 2e-16 *** --- Signif. codes: 0 *** ** 0.01 * Residual standard error: 5.74 on 44 degrees of freedom Multiple R-squared: , Adjusted R-squared: 0.95 F-statistic: on 5 and 44 DF, p-value: < 2.2e-16

18 ß

19 ß 0 0 X X ß

20 GA_recomb<-function(G) { npop=nrow(g); nvar=ncol(g); child = G; # parent 1 G.permuted = G[sample(1:nPop),] # parent 2 recomb.idx = which(matrix(sample(0:1,npop*nvar,replace=t), nrow=npop)==1) # crossover idx child[recomb.idx] = G.permuted[recomb.idx] return(child)

21 GA_mutate = function(child, p){ n.pop = nrow(child) n.var = ncol(child) mut.mat= matrix((runif(n.pop*n.var) < p), nrow = n.pop) child = abs(child-mut.mat) # abs(0-1) -> 1; abs(1-1) -> 0 # 0 -> 1; 1 -> 0 return(child)

22 GA_survive<-function(G, child, fitg, fitc){ npop=nrow(g); fitt=c(fitg,fitc); fitmax=sort(fitt, decreasing = TRUE, index.return=true) # maximize adjst. R^2 tempx=rbind(g,child); G=tempX[fitMax$ix[1:nPop],] return(g)

23 See courselog for MAIN function # RESULT > res<-ga(g,100,0.1,x,y) > head(res$g) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] [2,] [3,] [4,] [5,] [6,]

24 niv = GHI _K"L#EDMLKN" A KO N CDE i model x1 x2 x3 x1*x2 x1*x3 x2*x3 x1*x2*x

25 T niv = A 10 i KO = 637

26 Combn(Nvar, Nint) > combn(3,2) [,1] [,2] [,3] [1,] [2,] > combn(3,3) [,1] [1,] 1 [2,] 2 [3,] 3

27 C = combn(3,2) > paste("x",c,sep="") [1] "X1" "X2" "X1" "X3" "X2" "X3 > temp.lab =tostring(paste("x",c[, i.comb], sep = "")) [1] "X1, X2" > gsub(",", "*", temp.lab) [1] "X1* X2 n.comb = ncol(c) for (i.comb in 1: n.comb){ temp.label = tostring(paste("x",c[, i.comb], sep = "")) var.labels = c(var.labels, gsub(",", "*", temp.label) ) > var.labels [1] "X1* X2" "X1* X3" "X2* X3"

28 mk.labels <- function(c){ var.labels = c() n.comb = ncol(c) for (i.comb in 1: n.comb){ temp.label = tostring(paste("x",c[, i.comb], sep = "")) var.labels = c(var.labels, gsub(",", "*", temp.label) ) return(var.labels)

29 n.var = 3 # single variable var.labels = paste("x", 1:n.var, sep = "") # interaction terms for (i.interaction in 2:max.interaction ){ combination = combn(n.var, i.interaction) var.labels = c(var.labels, mk.labels(combination)) > var.labels [1] "X1 "X2" "X3" "X1* X2" "X1* X3" "X2* X3" "X1* X2* X3"

30 # creating formula for Linear Model model.def = paste("y ~", gsub(",", "+", tostring(var.labels[idx == 1]))) # example > model.def = paste("y ~", gsub(",", "+", tostring(var.labels[c(1,1,1,1,1,0,0) == 1]))) > model.def [1] "Y ~ X1+ X2+ X3+ X1* X2+ X1* X3 # running Linear Model model.lm = lm(model.def, data)

31 set.seed(20) nvar = 3 ndata = 50 X=matrix(rnorm(nVar*nData,mean=10,sd=5),ncol=nVar); y=x%*%c(1,1,-2)+rnorm(ndata,mean=0,sd=3); dat = data.frame(y=y,x1=x[,1],x2=x[,2],x3=x[,3]) # example model.def1 = paste("y ~", gsub(",", "+", tostring(var.labels[c(1,1,1,0,0,0,0) == 1]))) model.lm = lm(model.def1, dat)

32 # example model.def1 = paste("y ~", gsub(",", "+", tostring(var.labels[c(1,1,1,0,0,0,0) == 1]))) model.lm = lm(model.def1, dat) > summary(model.lm) Coefficients: Estimate Std. Error t value Pr(> t ) (Intercept) X e-15 *** X e-14 *** X < 2e-16 *** --- Signif. codes: 0 *** ** 0.01 * Residual standard error: on 46 degrees of freedom Multiple R-squared: , Adjusted R-squared: F-statistic: on 3 and 46 DF, p-value: < 2.2e-16

33 model IV1 IV2 IV3 IV4 IV5 IV6 IV7 IV8 IV9 IV10 1 c11 c12 c13 c14 c15 c16 c17 c18 c19 c110 2 c21 c22 c23 c24 c25 c26 c27 c28 c29 c210 Npop

34 set.seed(20); ndata = 100 X=matrix(rnorm(9*nData,mean=10,sd=2),ncol=9);X=cbind(rep(1,nData),X) y=x%*%c(10,2,5,-3,-5,0,0,0,0,0)+rnorm(ndata,mean=0,sd=2); > summary(lm(y~x[,2:10])) Coefficients: Estimate Std. Error t value Pr(> t ) (Intercept) ** X[, 2:10] < 2e-16 *** X[, 2:10] < 2e-16 *** X[, 2:10] < 2e-16 *** X[, 2:10] < 2e-16 *** X[, 2:10] X[, 2:10] X[, 2:10] X[, 2:10] X[, 2:10]

36 θ MVKWX K = Y θ 2DE#"L K if UNI K < 0.5 θ 2DE#"L] K otherwise σ K MVKWX = ] σ 2DE#"L 2DE#"L] K + σ K

37 σ K σ K exp τn 0,1 θ K θ K + N 0, σ K exp τn 0,1

38 θ MVKWX K = Y θ 2DE#"L K if UNI K < 0.5 θ 2DE#"L] K otherwise σ MVKWX K = 1 2 σ 2DE#"L 2DE#"L] K + σ K ES_recomb<-function(G) { nparent=nrow(g$x); nvar=ncol(g$x); child = G for (i_child in 1:nChild) { parentid=sample(1:nparent,2) coid=sample(c(0,1),nvar,replace=t) child$x[i_child,]=g$x[parentid[1],] child$x[i_child,which(coid==1)]=g$x[parentid[2],which(coid==1)] child$sigma[i_child,]=0.5*(g$sigma[parentid[1],]+g$sigma[parentid[2],]) return(child)

39 σ K σ K exp τn 0,1 θ K θ K + N 0, σ K ES_mutate<-function(child,tau){ nchild=nrow(child$x);nvar=ncol(child$x) child$sigma<-child$sigma*exp(matrix(rnorm(nchild*nvar)*tau,nrow=nchild)) child$x=child$x+child$sigma*matrix(rnorm(nchild*nvar),nrow=nchild,ncol=nvar) return(child)

40 ES_survive<-function(G, child, fitg, fitc){ nparent=nrow(g$x); fitt=c(fitg,fitc); fitmin=sort(fitt,index.return=t) tempx=rbind(g$x,child$x); temps=rbind(g$sigma,child$sigma) G$x=tempX[fitMin$ix[1:nParent],] G$sigma=tempS[fitMin$ix[1:nParent],] return(g)

42 set.seed(20); ndata = 100 X=matrix(rnorm(9*nData,mean=10,sd=2),ncol=9);X=cbind(rep(1,nData),X) y=x%*%c(10,2,5,-3,-5,0,0,0,0,0)+rnorm(ndata,mean=0,sd=2); fun04<-function(b,x,y){ yhat<-x%*%b return(sum((y-yhat)^2))

43 > res=es(g, fun04, 10000, 2,X,y) > head(round(res$g$x,3)) [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10] [1,] [2,] > summary(lm(y~x[,2:10])) Coefficients: Estimate Std. Error t value Pr(> t ) (Intercept) ** X[, 2:10] < 2e-16 *** X[, 2:10] < 2e-16 *** X[, 2:10] < 2e-16 *** X[, 2:10] < 2e-16 *** X[, 2:10] X[, 2:10] X[, 2:10] X[, 2:10] X[, 2:10]

45 v K t = v K t 1 + w K 2 U pb K θ K t 1 + w K l U gbk θ K t 1 θ K t θ K t 1 + v K t

47 # initialization swarm.hist = array(0, c(nrow(g), ncol(g), maxiter)) swarm.hist[,,1]=g p.b.hist = apply(g,1,func) global.best.v = min(p.b.hist) p.best = G g.best =matrix(g[which.min(p.b.hist),], nrow=nrow(g), ncol=ncol(g), byrow=t) v = matrix(0,nrow = nrow(g),ncol = ncol(g))

48 # main loop for (i.iter in 2:maxIter){ v = v + wp*runif(1)*(p.best - G) + wg*runif(1)*(g.best - G) G = G + v fitg = apply(g,1,func) if (min(fitg) < global.best.v){ g.best = matrix(g[which.min(fitg),],nrow=nrow(g),ncol=ncol(g),byrow=t) global.best.v = min(fitg) idx = which(fitg < p.b.hist) p.best[idx,] = G[idx,] p.b.hist= fitg swarm.hist[,,i.iter]=g

49 z = x n 16x ] 5x + y n 16y ] 5y

51 y = p 0 (w x + w ] x ] θ) 1 (w x + w ] x ] > θ)

52 y = p 0 (w x + w ] x ] θ) 1 (w x + w ] x ] > θ). y = p 0 (0.5x + 0.5x ] 0.7) 1 (0.5x + 0.5x ] > 0.7)

53 y = p 0 ( 0.5x 0.5x ] 0.7) 1 ( 0.5x 0.5x ] > 0.7) y = p 0 (0.5x + 0.5x ] 0.3) 1 (0.5x + 0.5x ] > 0.3)

54 AND.gate <- function(x1, x2){ w1 = 0.5 w2 = 0.5 theta = 0.7 y.temp = w1*x1 + w2*x2 if (y.temp <= theta){ y = 0 else { y = 1 return(y)

55 AND.gate <- function(x1, x2){ w1 = 0.5; w2 = 0.5; theta = 0.7 return(as.numeric(w1*x1 + w2*x2 > theta))

56 NAND.gate <- function(x1, x2){ w1 = -0.5; w2 = -0.5; theta = -0.7 return(as.numeric(w1*x1 + w2*x2 > theta)) OR.gate <- function(x1, x2){ w1 = 0.5; w2 = 0.5; theta = 0.3 return(as.numeric(w1*x1 + w2*x2 > theta))

57 y = p 0 (w x + w ] x ] θ) 1 (w x + w ] x ] > θ) y = p 0 (b + w x + w ] x ] 0) 1 (b + w x + w ] x ] > 0)

58 AND.gate <- function(x1, x2){ w1 = 0.5 w2 = 0.5 b = -0.7 y.temp = w1*x1 + w2*x2 + b if (y.temp <= 0){ y = 0 else { y = 1 return(y)

59 AND, NAND, and OR AND.gate <- function(x1, x2){ w1 = 0.5; w2 = 0.5; b = -0.7 return(as.numeric(w1*x1 + w2*x2 + b > 0)) NAND.gate <- function(x1, x2){ w1 = -0.5; w2 = -0.5; b = 0.7 return(as.numeric(w1*x1 + w2*x2 + b > 0)) OR.gate <- function(x1, x2){ w1 = 0.5; w2 = 0.5; b = -0.3 return(as.numeric(w1*x1 + w2*x2 + b > 0))

62 XOR.gate <- function(x1, x2){ gate1 <- NAND.gate(x1,x2) gate2 <- OR.gate(x1,x2) y <- AND.gate(gate1,gate2) return(y)

63 y = p 0 (b + w x + w ] x ] 0) 1 (b + w x + w ] x ] > 0) y = h b + w x + w ] x ] 0 (x 0) h x = p 1 (x > 0)

64 a = b + w x + w ] x ] y = h a

65 0 (x 0) h x = p 1 (x > 0) step.func <- function(x){ return(as.numeric(x > 0))

66 0 (x 0) h x = p 1 (x > 0) step.func <- function(x){ return(as.numeric(x > 0)) x = seq(-5, 5, 0.1) y = step.func(x) plot(x,y, ylab = 'y', xlab = 'a', type ="l", lwd =2)

67 h x = 0exp w1 sigmoid.func <- function(x){ return(1(1+exp(-x)))

68 h x = 0exp w1 sigmoid.func <- function(x){ return(1(1+exp(-x))) y = sigmoid.func(x) plot(x,y, ylab = 'y', xlab = 'a', type ="l", lwd =2)

69 y = sigmoid.func(x) y.step = step.func(x) y.sigm = sigmoid.func(x) plot(x,y.step, ylab = 'y', xlab = 'a', type ="l", lwd =2) lines(x,y.sigm, lwd =2, lty = 2)

70 0 (x 0) h x = p x (x > 0) relu.func <- function(x){ return(pmax(0,x)) y.relu = relu.func(x) plot(x,y.relu, ylab = 'y', xlab = 'a', type ="l", lwd =2)

71 a a ] a ] a ]] b b ] = a b + a ] b ] a b ] + a ] b ]] b ] b ]] a ] b + a ] b ] a ] b ] + a ]] b ]] = A = matrix(1:4, nrow = 2, byrow = T) B = matrix(5:8, nrow = 2, byrow = T) > A%*%B [,1] [,2] [1,] [2,] 43 50

72 A = matrix(1:6, nrow = 3, byrow = T) B = matrix(7:8, nrow = 2, byrow = T > A%*%B [,1] [1,] 23 [2,] 53 [3,] 83

73 b w w ] a a ] a z z ] z a ] a ] ] z ] z ] ] w 緑 : 層 : 次層のニューロン : 前層のニューロン

74 b w w ] a a ] a z z ] z X = x x ] A 1 = a a ] a B 1 = b b ] b W 1 = w w ] w w ] w ]] w ] A 1 = XW 1 + B Z 1 = z z ] z = h A 1 x = c(1,0.5) W1 = matrix((1:6)*0.1, nrow = 2) B1 = (1:3)*0.1 > A1 [,1] [,2] [,3] [1,] Z1 = sigmoid.func(a1) > Z1 [,1] [,2] [,3] [1,]

75 a z a ] z ] a ] a z ] z a ] ] z ] ] w 緑 : 層 : 次層のニューロン : 前層のニューロン

76 x = c(1,0.5) W1 = matrix((1:6)*0.1, nrow = 2) B1 = (1:3)*0.1 Z1 = sigmoid.func(a1) W2 = matrix((1:6)*0.1, nrow = 3) B2 = c(0.1, 0.2) A2 = Z1%*%W2 + B2 Z2 = sigmoid.func(a2) W3 = matrix((1:4)*0.1, nrow = 2) B3 = c(0.1, 0.2) A3 = Z2%*%W3+ B3 Z3 = A3 > Z3 [,1] [,2] [1,] X = x x ] A 1 = a a ] a B 1 = b b ] b W 1 = w w ] w w ] w ]] w ] A 1 = XW 1 + B 1 Z 1 = h A 1 A 2 = Z 1 W 2 + B 2 Z 2 = h A 2 A 3 = Z 2 W 3 + B 3 Z 3 = h A 3

77 # function to initialize 3L network init.3l.network <- function(){ W1 = matrix((1:6)*0.1, nrow = 2) B1 = (1:3)*0.1 W2 = matrix((1:6)*0.1, nrow = 3) B2 = c(0.1, 0.2) W3 = matrix((1:4)*0.1, nrow = 2) B3 = c(0.1, 0.2) return(list(w1 = W1, B1 = B1, W2 = W2, B2 = B2, W3 = W3, B3 = B3)) # feedforward process forward.3l <- function(network, x){ A1 = x%*%network$w1 + network$b1 Z1 = sigmoid.func(a1) A2 = Z1%*%network$W2 + network$b2 Z2 = sigmoid.func(a2) A3 = Z2%*%network$W3 + network$b3 Z3 = sigmoid.func(a3) A3 = Z3 return(a3) network<-init.3l.network() y = forward.3l(network, c(1, 0.5)) > y [,1] [,2] [1,]

78 σ( ) σ( ) y = a y = exp D ƒ exp D σ( ) y = exp D ƒ0 ˆ exp D 0 ˆ σ( )

79 y = exp D ƒ0 ˆ exp D 0 ˆ > exp(a)sum(exp(a)) [1] NaN NaN NaN σ( ) softmax.func <- function(x){ max.x = max(x) return(exp(x-max.x)sum(exp(x-max.x))) > softmax.func(a) [1] e e e-09

80 train <- read.csv(' header=true) train <- data.matrix(train) train.x <- train[,-1] train.y <- train[,1] train.x <- t(train.x255) download.file(" "trnetwork.rdata") load("trnetwork.rdata") network=trnetwork

81 n.train = ncol(train.x) correct.cl = 0 conf.matrix = matrix(0,10,10) for (i.loop in 1:n.train){ y = forward.3l(network,train.x[,i.loop]) max.y = max.col(y) conf.matrix[max.y, (train.y[i.loop]+1)] = conf.matrix[max.y, (train.y[i.loop]+1)] + 1 accuracy = sum(diag(conf.matrix))n.train

83 batch_size = 200 conf.matrix = matrix(0,10,10) for (i.batch in seq(1,n.train, batch_size)){ y = forward.3l(network, train.x[,(i.batch:(i.batch+batch_size-1))]) pred = max.col(y) conf.matrix = conf.matrix + table(pred,(train.y[i.batch:(i.batch+batch_size-1)]+1)) accuracy = sum(diag(conf.matrix))n.train

84 Batch > system.time({ + for (i.batch in seq(1,n.train, batch_size)){ + y = forward.3l(network, train.x[,(i.batch:(i.batch+batch_size-1))]) + pred = max.col(y) + conf.matrix = conf.matrix+table(pred,(train.y[i.batch:(i.batch+batch_size-1)]+1)) + + ) user system elapsed Online > system.time({ + for (i.loop in 1:n.train){ + y = forward.3l(network,train.x[,i.loop]) + max.y = max.col(y) + conf.matrix[max.y, (train.y[i.loop]+1)] = conf.matrix[max.y, (train.y[i.loop]+1)] ) user system elapsed

DAA09

DAA09 > summary(dat.lm1) Call: lm(formula = sales ~ price, data = dat) Residuals: Min 1Q Median 3Q Max -55.719-19.270 4.212 16.143 73.454 Coefficients: Estimate Std. Error t value Pr(> t ) (Intercept) 237.1326