Posted on May 16, 2017 by matsukat

データ解析基礎論A Week05 t検定

dat<-read.csv("http://www.matsuka.info/data_folder/datWA01.txt")
mean.M <-mean(dat$h[dat$gender=="M"])
sigma = 10
n.M = length(dat$h[dat$gender=="M"])
z.value=(mean.M-171)/(sqrt(sigma/n.M))
(1-pnorm(abs(z.value)))*2

ssize = c(24,25,26,23.5,25,27,24,22,27.5,28)
ssize.mean = mean(ssize)
ssize.var = var(ssize)
N = 10
t.value=(ssize.mean-24)/(sqrt(ssize.var/N))
(1-pt(abs(t.value),df=9))*2

h.mean.M <-mean(dat$h[dat$gender=="M"])
h.var.M <- var(dat$h[dat$gender=="M"])
n.M = length(dat$h[dat$gender=="M"])
t.value=(h.mean.M-171)/(sqrt(h.var.M/n.M))
(1-pt(abs(t.value),df = (n.M-1)))*2

# RCMD
A=c(12,19,10,10,14,18,15,11,16)
B=c(15,20,16,14,17,16,12,12,19)
d=A-B
tValue<-mean(d)/sqrt(var(d)/length(d))
(1-pt(abs(tValue), df=8))*2

t.test(d,mu=0)


X1=c(78,70,66,76,78,76,88,76)
X2=c(76,72,60,72,70,72,84,70)
t.value=(mean(X1)-mean(X2))/sqrt((var(X1)+var(X2))/8)
2*(1-pt(abs(t.value),14))

Posted on May 10, 2017 by matsukat

認知情報解析　ch04

MSE <- function(target, y){
  return(0.5*sum((target-y)^2))
}

t = rep(0,10)
t[3]=1
y = c(0.1, 0.05, 0.6, 0, 0.05, 0.1, 0, 0.1, 0, 0)

x = seq(0,1,0.01)
plot(x,-log(x),lwd = 2)
cross.entropy = function(y, target){
  delta = 1e-7;
　R = nrow(as.matrix(y))
  return(-sum(target*log(y + delta))/R)
}

numerical.diff = function(func, x){
  h = 1e-4
  plusH = do.call(func,list(x+h))
  minusH = do.call(func,list(x-h))
  num.diff = (plusH - minusH)/(2*h)
  return(num.diff)
}

func01 = function(x){
  return(0.01*x^2+0.1*x)
}

x = seq(0,20,0.1)
y = func01(x)
plot(x,y,xlab ="x", ylab = "f(x)",type = "l",lwd =2)
ND.5 = numerical.diff('func01',5)
abline(a = func01(5)-ND.5*5, b = ND.5, col = 'red',  lwd =2)
abline(v = 5, lty = 2, col = 'red')
ND.10 = numerical.diff('func01',10)
abline(a = func01(10)-ND.10*10, b = ND.10, col = 'blue',lwd = 2)
abline(v = 10, lty = 2, col = 'blue')


func02 = function(x0, x1){
  return(x0^2 + x1^2)
}
func02.x0 = function(x0){
  return(x0^2)
}
func02.x1 = function(x1){
  return(x1^2)
}

func02R = function(x){
  return(x[1]^2 + x[2]^2)
}
numerical.grad <- function(func, x){
  h = 1e-4
  R = nrow(x)
  C = ncol(x)
  grad = matrix(0, R, C)
  for (i.col in 1:C){
    for (i.row in 1:R){
      temp.x = x[i.row,i.col]
      x[i.row, i.col] = temp.x + h
      plusH = do.call(func, list(x))
      x[i.row, i.col] = temp.x - h
      minusH = do.call(func,list(x))
      grad[i.row, i.col] = (plusH - minusH)/(2*h)
      x[i.row, i.col] = temp.x
    }
  }
  return(grad)
}

numerical.grad("func02R",matrix(c(3,4),nrow=1))
numerical.grad("func02R",matrix(c(0,4),nrow=1))
numerical.grad("func02R",matrix(c(3,0),nrow=1))

require(plot3D)
x = seq(-2,2,0.2)
y = seq(-2,2,0.2)
M = mesh(x,y)
R = nrow(M$x)
C = nrow(M$x)
scaling = 0.05
plot(c(),c(),xlim = c(-2,2),ylim=c(-2,2))
for (i.col in 1:C){
  for (i.row in 1:R){
    ng = numerical.grad("func02R",matrix(c(M$x[i.row,i.col],M$y[i.row,i.col]),nrow=1))
    arrows(M$x[i.row,i.col],M$y[i.row,i.col],
           (M$x[i.row,i.col]-ng[1]*scaling),(M$y[i.row,i.col]-ng[2]*scaling),
           length = 0.05)
  }
}

grad.desc <- function(func, init.x, lr, n.iter){
  x = init.x
  for (i.iter in 1:n.iter) {
    grad = numerical.grad(func, x)
    x = x - lr*grad
  }
  return(x)
}
x.init = matrix(c(-3,4),nrow = 1)
grad.desc("func02R",x.init,0.1,100)

x = seq(-4,4,0.2)
y = seq(-4,4,0.2)
M = mesh(x,y)
Z = as.vector(M$x^2)+as.vector(M$y^2) 
Z.mesh = matrix(Z,nrow(M$x))
contour(x,y,Z.mesh,drawlabels = F)
grad.desc2 <- function(func, init.x, lr, n.iter){
  x = init.x
  x.hist = init.x
  for (i.iter in 1:n.iter) {
    grad = numerical.grad(func, x)
    x = x - lr*grad
    x.hist = rbind(x.hist,x)
  }
  return(x.hist)
}
gd = grad.desc2("func02R",x.init,0.1,100)
points(gd,col = 'green',pch=20)

# manual implementation
w = matrix(c(0.47355232,0.85557411,0.9977393,0.03563661,0.84668094,0.69422093),nrow=2)
x = matrix(c(0.6, 0.9), nrow=1)
t = c(0,0,1)
nn.predict <- function(w,x){
  return(x%*%w)
}

loss.func = function(w, x, t){
  pred = nn.predict(w,x)
  y = softmax.func(pred)
  return(cross.entropy(y, t))
}


numerical.gradCE <- function(func, w, x, t){
# input args
# func: name of function
# w   : weight
# x   : input 
# t   : target output
##############################################
  h = 1e-4
  R = nrow(w)
  C = ncol(w)
  grad = matrix(0, R, C)
  for (i.col in 1:C){
    for (i.row in 1:R){
      temp.w = w[i.row,i.col]
      w[i.row, i.col] = temp.w + h
      plusH = do.call(func, list(w,x,t))
      w[i.row, i.col] = temp.w - h
      minusH = do.call(func,list(w,x,t))
      grad[i.row, i.col] = (plusH - minusH)/(2*h)
      w[i.row, i.col] = temp.w
    }
  }
  return(grad)
}

dW = numerical.gradCE("loss.func",w,x,t)

### ch 4.5 2-layer NN  ###
init.2LN <- function(n.input, n.hidden, n.output, w.std = 0.01){
  W1 = matrix(rnorm(n.input*n.hidden,0,w.std),nrow = n.input)
  B1 = matrix(rnorm(n.hidden,0,w.std),nrow =1)
  W2 = matrix(rnorm(n.hidden*n.output,0,w.std),nrow = n.hidden)
  B2 = matrix(rnorm(n.output,0,w.std),nrow =1)
  return(list(W1 = W1, B1 = B1, W2 = W2, B2 = B2))
}

softmax.2LN <- function(x){
  max.x = apply(x,1,max)
  C = ncol(x)
  x = x - max.x%*%matrix(1,nrow=1,ncol=C)
  return(exp(x)/rowSums(exp(x)))
}

sigmoid.func <- function(x){
  return(1/(1+exp(-x)))
}

pred.2LN <- function(params, x){
  NR = nrow(x)
  a1 = x%*%params$W1 + matrix(1,nrow = NR)%*%params$B1
  z1 = sigmoid.func(a1)
  a2 = z1%*%params$W2 + matrix(1,nrow = NR)%*%params$B2
  y = softmax.2LN(a2)
  return(y)
}

loss.2LN = function(params, x, t){
  y = pred.2LN(params,x)
  return(cross.entropy(y, t))
}

numerical.grad2LN <- function(func, params, x, t) {
  # input args
  # func: name of function
  # w   : weight
  # x   : input
  # t   : target output
  ##############################################
  h = 1e-4; n.list = length(params); grad = params
  for (i.list in 1:n.list) {
    R = nrow(params[[i.list]])
    C = ncol(params[[i.list]])
    grad[[i.list]] = matrix(0, R, C)
    for (i.col in 1:C) {
      for (i.row in 1:R) {
        temp.w = params[[i.list]][i.row, i.col]
        params[[i.list]][i.row, i.col] = temp.w + h
        plusH = do.call(func, list(params, x, t))
        params[[i.list]][i.row, i.col] = temp.w - h
        minusH = do.call(func, list(params, x, t))
        grad[[i.list]][i.row, i.col] = (plusH - minusH) / (2 * h)
        params[[i.list]][i.row, i.col] = temp.w
      }
    }
  }
  return(grad)
}

## example using IRIS data set
train.x = as.matrix(iris[,1:4])
train.y.temp = as.numeric(iris[,5])
train.y = matrix(0,nrow = nrow(train.x), ncol =3)
train.y[which(train.y.temp==1), 1]=1
train.y[which(train.y.temp==2), 2]=1
train.y[which(train.y.temp==3), 3]=1

params = init.2LN(4,15,3,0.01)
batch_size = 7; n.iter =2000; lambda =0.05
n.train = nrow(train.x)
loss = rep(0,n.iter)
for (i.iter in 1:n.iter){
  batch_mask = sample(1:n.train, batch_size)
  x.batch = train.x[batch_mask,]
  t.batch = train.y[batch_mask,]
  dW = numerical.grad2LN("loss.2LN",params,x.batch,t.batch)
  params$W1 = params$W1 - lambda*dW$W1
  params$B1 = params$B1 - lambda*dW$B1
  params$W2 = params$W2 - lambda*dW$W2
  params$B2 = params$B2 - lambda*dW$B2
  loss[i.iter] = loss.2LN(params,x.batch,t.batch)
}

Posted on May 10, 2017 by matsukat

認知情報解析・課題１・解答例

relu.func <- function(x){
 y = apply(x,2,function(x) pmax(x,0))
 return(y)
}

sigmoid.func <- function(x){
  return(1/(1+exp(-x)))
}

init.network <- function(n.neurons){
  n.layer = length(n.neurons)
  W = list(); b = list()
  for (i.layer in 1:(n.layer-1)){
    W[[i.layer]] = 
      matrix(rnorm(n.neurons[i.layer]*n.neurons[(i.layer+1)]),nrow=n.neurons[i.layer])
    b[[i.layer]] =  matrix(rnorm(n.neurons[(i.layer+1)]), nrow = 1)
  }
  return(list(W = W,b = b))
}

activation <- function(A, actFun){
  if (actFun == "sigmoid"){
    return(sigmoid.func(A))
  }
  if (actFun == "relu"){
    return(relu.func(A))
  }
}

feedforward <- function(network, x, actFun) {
  n.layer <- length(network$W)
  batch.size = nrow(x)
  for (i.layer in 1:n.layer){
    A = x%*%network$W[[i.layer]] + 
      matrix(1,nrow=batch.size,ncol = 1)%*%network$b[[i.layer]]
    x = activation(A, actFun[i.layer])
  }
  return(x)
}

Posted on May 8, 2017 by matsukat

データ解析基礎論a Week04

nSample=10;nRep=10^5;
CLT.unif <- function(nSample, nRep) {
  x=matrix(runif(nSample*nRep),nrow=nSample);
  x.means<-colMeans(x)
  hist(x.means,50,main='Distribution of Means of Uniform Distribution', 
    xlab='Means', probability=T)
  x.means.dens<-density(x.means)
  lines(x.means.dens,lwd=3,lty=1,col='red')
  x2=seq(0,1,0.001);CLT=dnorm(x2,mean=0.5,sd=(sqrt(1/12))/(sqrt(nSample)))
  lines(x2,CLT,lwd=3,lty=3,col='cyan')
  legend("topright",c("Density of Actual Means","Normal Distribution"), 
   lty=c(1,3), col=c('red','cyan'),lwd=3)
}

> CLT.unif(10,100000)

Posted on April 25, 2017 by matsukat

データ解析基礎論A Week03

dat<-read.table("http://www.matsuka.info/data_folder/aov01.txt",header=T)
head(dat)
summary(dat)
boxplot(dat$shoesize, col = "skyblue", main = "Dist. of Shoesize",
 ylab = "Size of shoe (in cm)")
boxplot(dat$h, col = "coral", main = "Dist. of Height", ylab = "Height (in cm)")

dat.meter = dat[,1:2]*0.01 
dat.h.ext5 = dat$h+5
dat.ss.ext1 = dat$shoesize+1
hANDshoe = dat$h+dat$shoesize
dat.h.meter.ext5 = dat$h*0.01+0.05

Posted on April 19, 2017 by matsukat

ch03 neural networks

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)

sigmoid.func <- function(x){
  return(1/(1+exp(-x)))
}

y = sigmoid.func(x)
plot(x,y, ylab = 'y', xlab = 'a', type ="l", lwd =2)

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)

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)

A = matrix(1:4, nrow = 2, byrow = T)
B = matrix(5:8, nrow = 2, byrow = T)

A = matrix(1:6, nrow = 3, byrow = T)
B = matrix(7:8, nrow = 2, byrow = T)

x = c(1,0.5)
W1 = matrix((1:6)*0.1, nrow = 2)
B1 = (1:3)*0.1
A1 = x%*%W1 + B1
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

# 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))

a = c(1010,1000,990)
exp(a)/sum(exp(a))

softmax.func <- function(x){
  max.x = max(x)
  return(exp(x-max.x)/sum(exp(x-max.x)))
}
  

train <- read.csv('http://peach.l.chiba-u.ac.jp/course_folder/MNSTtrain.csv', 
  header=TRUE)
train <- data.matrix(train)
train.x <- train[,-1]
train.y <- train[,1]
train.x <- t(train.x/255)
download.file("http://peach.l.chiba-u.ac.jp/course_folder/trNetwork.Rdata",
  "trNetwork.Rdata")
load("trNetwork.Rdata")
network=trNetwork

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

# batch 
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

Posted on April 17, 2017 by matsukat

データ解析基礎論A グラフの基礎

データ解析基礎論A 第2週　グラフの基礎

# descriptive statistics
dat<-read.table("http://www.matsuka.info/data_folder/aov01.txt",header=T)
summary(dat)
mean(dat$shoesize)
var(dat[,1:2])
cov(dat[,1:2])
cor(dat[,1:2])

# basics - plotting
x=seq(-3,3,0.1);y=x^2;
plot(x,y)
plot(x,y,col='red')
plot(x,y,pch=20)
plot(x,y,type='l')
plot(x,y,type='l',lty=4,lwd=3)
plot(x,y,main="THIS IS THE TITLE", xlab="Label for X-axis",ylab="Label for Y-axis")
plot(x,y,main="TITLE", xlab="X here",ylab="Y here",xlim=c(-3.5,3.5),ylim=c(-0.5, 10))
plot(x,y,col='blue',type='o',lty=2,pch=19,lwd=3,main="Y=X*X", xlab="X",ylab="X*X",
  xlim=c(-3.5,3.5),ylim=c(-0.5, 10))

# histogram
dat<-read.table("http://www.matsuka.info/data_folder/aov01.txt",header=T)
hist(dat$h,main="Histogram of Height",xlab="Height",col='blue',xlim=c(140,190))
dens<-density(dat$h);
hist(dat$h,main="Histogram of Height",xlab="Height",xlim=c(140,190),probability=T)
lines(dens,lwd=2,col='red',lty=2)

# boxplot
boxplot(dat$h,main="Boxplot of Height",ylab="Height",col='cyan',ylim=c(140,190))
boxplot(dat$h~dat$gender,main="Distribution of Height by Gender",ylab="Gender",
 xlab="Height",col=c('blue','cyan'),ylim=c(140,190),horizontal=T)
boxplot(dat$h~dat$gender+dat$affil, main="Distribution of Height by Gender and Affiliation",
  ylab="Gender x Affiliation", xlab="Height", col=c('blue','cyan','red','magenta'),
  ylim=c(140,190),horizontal=T)

# barplot
install.packages("gplots")
library(gplots)

means <- tapply(dat$h, dat$gender, mean)
sds<-tapply(dat$h,dat$gender,sd)     
ns<-tapply(dat$h,dat$gender,length) 
sems = sds/sqrt(ns) 
barplot2(means, plot.ci=T,
 ci.l = means - sems,
 ci.u = means + sems,
 ylim = c(140,180),
 names.arg = c("Female","Male"),
 xpd = F,
 ylab = "height",
 xlab = "gender")

# another barplot 
means <- tapply(dat$h,list(dat$gender,dat$affil),mean)
sds <- tapply(dat$h,list(dat$gender,dat$affil),sd)
ns <- tapply(dat$h,list(dat$gender,dat$affil),length)
sem = sds/sqrt(ns)
barplot2(means[1:4], plot.ci=T, ci.l=means[1:4]-sem[1:4],   
        ci.u=means[1:4] + sem[1:4], ylim=c(150,175), 
    names.arg=c("Female,CS","Male,CS","Female,PSY","Male,PSY"),   
        xpd=F,ylab="height",xlab="gender & affiliation")

# histogram again
par(mfrow=c(1,2))
hist(dat[dat$gender=='F',]$h,main="Dist. of Height for Female Participants",
 xlab="Height",xlim=c(140,190),probability=T)
dens.F=density(dat[dat$gender=='F',]$h);lines(dens.F,col='blue',lwd=2)
hist(dat[dat$gender=='M',]$h,main="Dist. of Height for Male Participants",
 xlab="Height",xlim=c(140,190),probability=T,ylim=c(0,0.08))
dens.M=density(dat[dat$gender=='M',]$h);lines(dens.M,col='green',lwd=2)

par(mfrow=c(2,1))

par(mfrow=c(1,1))
plot(dens.F,col='blue',lwd=2,main="Dist. of Height by gender",xlab='Height',
  ylab='density',xlim=c(140,190))
lines(dens.M,col='green',lwd=2)
legend("topleft", c('Female','Male'),col=c('blue','green'),cex=1.5,lwd=2)

# inserting text
text(157.5,0.04,'Female',col='blue',cex=2)
text(170,0.04,'Male',col='green',cex=2)

# scatterplot
plot(dat$shoesize,dat$h,main="Relationship b/w shoesize and height",xlab='shoe size', 
  ylab='height',pch=19,col='red')
text(22,175,paste("r =",substr(cor(dat$shoesize,dat$h),1,5)),cex=1.5)
abline(h=mean(dat$h),col='blue');
abline(v=mean(dat$shoesize),col='green');
text(21.5,165,'mean height',col='blue')
text(25.7,145,'mean shoesize',col='green')
abline(lm(dat$h~dat$shoesize),lty=2,lwd=2)

plot(dat[dat$gender=='F',]$shoesize,dat[dat$gender=='F',]$h,
 main="Relationship b/w shoesize and height",xlab='shoesize',ylab='height',
 cex.lab=1.5,pch=19,col='blue',xlim=c(20,29),ylim=c(140,190))
lines(dat[dat$gender=='M',]$shoesize,dat[dat$gender=='M',]$h,type='p',pch=15,col='green')
legend("topleft", c('Female','Male'), pch =c(19,15),col=c('blue','green'),cex=1.5)
dat.reg<-read.csv("http://www.matsuka.info/data_folder/tdkReg01.csv", header=T)
plot(dat.reg,pch=20,col=c('blue'))
dat.pca<-read.table("http://www.matsuka.info/data_folder/tdkPCA01.txt",header=T)

# intro to central limit theorem
ckCLT=function(n_rep,n_sample){
  dat<-matrix(rnorm(n_rep*n_sample),nrow=n_rep,ncol=n_sample);
  means<-rowMeans(dat)}
n_rep=10^6
n5=ckCLT(n_rep,5)
hist(n5,main="Dist. of sample meanx",xlab="sample mean",xlim=c(-3,3),probability=T)
den5=density(n5);lines(den5,col='blue',lwd=2)

n10=ckCLT(n_rep,10)
n25=ckCLT(n_rep,25)
n100=ckCLT(n_rep,100)
plot(den5,col='blue',lwd=2,,main="Dist. of sample meanx",xlab="sample mean",
  xlim=c(-2,2),ylim=c(0,4))
den10=density(n10);lines(den10,col='red',lwd=2)
den25=density(n25);lines(den25,col='black',lwd=2)
den100=density(n100);lines(den100,col='green',lwd=2)
legend("topleft", c('N=5','N=10','N=25','N=100'),col=c('blue','red','black','green'),
  cex=1.5,lwd=2)

Posted on April 12, 2017 by matsukat

データ解析基礎論a WAA01

データ解析基礎論A　Weekly Assignment WAA01
提出期限：2017.04.18 4限開始まで
提出物：　Rのコマンドとその出力

WAA01.1
下記のデータをdata.frameのコマンドを用いてRに入力してください。

> dat
 id gender score
  1      M    73
  2      M    64
  3      M    78
  4      F    74
  5      F    84
  6      F    78

WAA01.2
gender == “M” や　gender == “F” などを用いて、上記のデータ男性と女性のscoreの平均値を求めてください。

WAA01.3
このファイルをダウンロードして、scoreが100の男性と女性の数を調べてください。

解答例
WAA01.1
例１
dat<-data.frame(id = 1:6, gender = c(rep("M",3),rep("F",3)), 
  score = c(73,64,78,74,84,78))
例２
id = 1:6
gender = c(rep("M",3),rep("F",3))
score = c(73,64,78,74,84,78)
dat2 <- data.frame(id = id, gender = gender, score = score)

WAA01.2
例１
mean(dat$score[which(dat$gender=="M")])
mean(dat$score[which(dat$gender=="F")])
例２
dat.male = dat[which(dat$gender=="M"),]
mean(dat.male$score)
dat.female = dat[which(dat$gender=="F"),]
mean(dat.female$score)

WAA01.3
例１
dat3 <- read.csv("http://peach.l.chiba-u.ac.jp/course_folder/waa01.csv")
head(dat3)
which(dat3[which(dat3$gender == "M"),]$score == 100)
which(dat3[which(dat3$gender == "F"),]$score == 100)
例２
dat3.male = dat3[which(dat3$gender == "M"), ]
dat3.female = dat3[which(dat3$gender == "F"), ]
dat3.male$id[which(dat3.male$score == 100)]
dat3.female$id[which(dat3.female$score == 100)]

Posted on April 12, 2017 by matsukat

認知情報解析 2017 – Ch 02

Chapter 2: Perceptron

# with THRESHOLD (theta)
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)
}

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

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))
}

# with BIAS (b)
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)
}

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))
}

NOR.gate <- function(x1, x2){
  w1 = -0.5; w2 = -0.5; b = 0.3
  return(as.numeric(w1*x1 + w2*x2 + b > 0))
}

plot.logic <- function(logic.oper){
  x1 = c(0,0,1,1);
  x2 = c(0,1,0,1);
  if (logic.oper == "and") {
    w1 = 0.5; w2 = 0.5; theta = 0.7;
    true.point = AND.gate(x1,x2)
  } else if (logic.oper == "or") {
    w1 = 0.5; w2 = 0.5; theta = 0.3;
    true.point = OR.gate(x1,x2)
  } else if (logic.oper == "nand") {
    w1 = -0.5; w2 = -0.5; theta = -0.7;
    true.point = NAND.gate(x1,x2)
  } else if (logic.oper == "nor"){
    w1 = -0.5; w2 = -0.5; theta = -0.3;
    true.point = NOR.gate(x1,x2)
  } else {warning("incompatible operator");stop() }
  plot(c(0,0,1,1),c(0,1,0,1),xlim = c(-0.5, 1.5), ylim = c(-0.5, 1.5), 
       pch = 20, cex= 2, col = true.point+1)
  abline(a = theta/w1, b = -w1/w2, lwd = 3)
}   
　
XOR.gate <- function(x1, x2){
  gate1 <- NAND.gate(x1,x2)
  gate2 <- OR.gate(x1,x2)
  y <- AND.gate(gate1,gate2)
  return(y)
}

plot.XOR <- function(){
  x1 = c(0,0,1,1);
  x2 = c(0,1,0,1);
  w11 = -0.5; w21 = -0.5; theta1 = -0.7
  w12 = 0.5; w22 = 0.5; theta2 = 0.3
  true.point = XOR.gate(x1, x2)
  plot(c(0,0,1,1),c(0,1,0,1),xlim = c(-0.5, 1.5), ylim = c(-0.5, 1.5), 
       pch = 20, cex= 2, col = true.point+1)
  abline(a = theta1/w11, b = -w11/w21, lwd = 3)
  abline(a = theta2/w12, b = -w12/w22, lwd = 3)    
}

Posted on April 11, 2017 by matsukat

データ解析基礎論a 2017 W01

データ解析基礎論a 2017 W01で使用するコードなど

data01<-data.frame(score = c(2,4,3,4),     
                   dose = c(rep(10,2),rep(100,2)),  
                   condition = rep(c('exp','control'),2))

dat01<-read.csv("http://www.matsuka.info/data_folder/temp_data01.txt",
                header=T)

dat02<-read.csv("http://www.matsuka.info/data_folder/temp_data02.txt", 
                header=T, row.name=1)

dat03<-read.table("http://www.matsuka.info/data_folder/temp_data03.txt",
                  header=T, row.name=4)

dat04<-data.frame(score=c(dat03$x,dat03$y,dat03$z),   
                  names=rep(c('sazae','masuo','tarachan'),3),
                  condition=sort(rep(c("x","y","z"),3)))

dat<-read.csv("http://www.matsuka.info/data_folder/datWA01.txt",   
              header=T);

mean(dat$shoesize[dat$gender == "M"])
mean(dat$shoesize[dat$gender == "F"])
mean(dat$shoesize[dat$h > 180])

CourseLog

Dept Cognitive & Info Science

Author Archives: matsukat