# 広域システム特別講義II 認知・社会モデル

```# cognitive modeling
ALC.init<-function(size) {
# size[1] = input dimension
# size[2] = number of exemplars
# size[3] = number of categories
alpha=matrix(1,nrow=1,ncol=size[1])/size[1]
w=matrix(rnorm(size[3]*size[2])*0.1,ncol=size[2])
c=2            # gain
phi=3          # decision strength
lrA=0.2        # learning rate for attentions
lrW=0.1        # learning rate for associations
return(list(alpha = alpha,
w = w,
lrA = lrA,
lrW = lrW,
c = c,
phi = phi))
}

# alcove forward process
alcove.forward <- function(parmSet, inp, exemplar){
diff= matrix(abs(matrix(1,n.exemp,1)%*%inp-exemplar),nrow=nrow(exemplar))
exemp=exp(-parmSet\$c*(diff%*%t(parmSet\$alpha)))
out=parmSet\$w%*%exemp
return(list(diff = diff, exemp = exemp, out = out))
}

# alcove backward process
alcove.backward <- function(res.forward,parmSet, target){
err=target-res.forward\$out
dW=parmSet\$lrW*res.forward\$exemp%*%t(err)
dA=-parmSet\$lrA*(t(err)%*%parmSet\$w*t(res.forward\$exemp))%*%res.forward\$diff*parmSet\$c
return(list(dW = dW, dA = dA))
}

### ALCOVE full implementation
ALC.init<-function(size) {
# size[1] = input dimension
# size[2] = number of exemplars
# size[3] = number of categories
alpha=matrix(1,nrow=1,ncol=size[1])/size[1]
w=matrix(rnorm(size[3]*size[2])*0.1,ncol=size[2])
c=2            # gain
phi=3          # decision strength
lrA=0.2        # learning rate for attentions
lrW=0.1        # learning rate for associations
return(list(alpha = alpha,
w = w,
lrA = lrA,
lrW = lrW,
c = c,
phi = phi))
}

# alcove forward process
alcove.forward <- function(parmSet, inp, exemplar){
diff= matrix(abs(matrix(1,n.exemp,1)%*%inp-exemplar),nrow=nrow(exemplar))
exemp=exp(-parmSet\$c*(diff%*%t(parmSet\$alpha)))
out=parmSet\$w%*%exemp
return(list(diff = diff, exemp = exemp, out = out))
}

# alcove backward process
alcove.backward <- function(res.forward,parmSet, target){
err=target-res.forward\$out
dW=parmSet\$lrW*res.forward\$exemp%*%t(err)
dA=-parmSet\$lrA*(t(err)%*%parmSet\$w*t(res.forward\$exemp))%*%res.forward\$diff*parmSet\$c
return(list(dW = dW, dA = dA))
}

# main function
alcove<-function(parmSet, inp, exemplar, targ,iter) {
# ----ALCOVE for R ---
# input arguments
#  parmSet - parameter set
# inp - input matrix (n_sample x n_dimension)
# exemplar - exemplar matrix (n_sample x n_dimension)
# targ - target matrix (n_sample x n_category)
# iter - # of training epochs
# ----ALCOVE for R ---
# initialization
inpDim = ncol(inp)
n.inp = nrow(inp)
n.exemp = nrow(exemplar)
n.cat = ncol(targ)
accu=rep(0,nrow=iter+1)
attn=matrix(0,nrow=iter+1,ncol=inpDim)
attn[1,]=alpha
# end initialization

# main loop
for (i_iter in 1:iter) {
ro=sample(1:n.inp, n.inp)
prob_temp=0;
for (i_tr in 1:n.inp) {
res.forward <- alcove.forward(parmSet, inp[ro[i_tr],], exemplar)
res.backward <- alcove.backward(res.forward, parmSet, targ[ro[i_tr],])
parmSet\$w = parmSet\$w+t(res.backward\$dW)
parmSet\$alpha = parmSet\$alpha+res.backward\$dA
parmSet\$alpha[which(parmSet\$alpha<0)]=0
pT=(exp(parmSet\$phi*res.forward\$out)/sum(exp(parmSet\$phi*res.forward\$out)))*targ[ro[i_tr],]
prob_temp=prob_temp+pT[which(!pT==0)]
}
accu[i_iter+1]=prob_temp/n.inp
attn[i_iter+1,]=alpha
}
attnN=attn/apply(attn,1, sum)
out=matrix(0,nrow=n.cat,ncol=n.inp)
for (i_tr in 1:n.inp) {
res.forward<-alcove.forward(parmSet, inp[i_tr,], exemplar)
out[,i_tr]=res.forward\$out
}
return(list(attn=attn,attnN=attnN,accu=accu,out=out,ps=parmSet))
}

exemplar = matrix(c(0,0,0,
0,0,1,
0,1,0,
0,1,1,
1,0,0,
1,0,1,
1,1,0,
1,1,1),byrow=T,nrow=8)
inp = exemplar

target1 = matrix(c(0,0,0,0,1,1,1,1,
1,1,1,1,0,0,0,0),nrow=8)

target2 = matrix(c(1,1,0,0,0,0,1,1,
0,0,1,1,1,1,0,0),nrow=8)

target3 = matrix(c(1,1,1,0,0,1,0,0,
0,0,0,1,1,0,1,1),nrow=8)

target4 = matrix(c(1,1,1,0,1,0,0,0,
0,0,0,1,0,1,1,1),nrow=8)

target5 = matrix(c(1,1,1,0,0,0,0,1,
0,0,0,1,1,1,1,0),nrow=8)

target6 = matrix(c(1,0,0,1,0,1,1,0,
0,1,1,0,1,0,0,1),nrow=8)

seed.num = 1;n.train = 50
set.seed(seed.num)
parmSet<-ALC.init(c(3,8,2))
result1<-alcove(parmSet,inp,exemplar,target1,n.train)
plot(result1\$accu,type='o',lwd=2,pch=20,col=1)
set.seed(seed.num)
parmSet<-ALC.init(c(3,8,2))
result2<-alcove(parmSet,inp,exemplar,target2,n.train)
lines(result2\$accu,type='o',lwd=2,pch=20,col=2)
set.seed(seed.num)
parmSet<-ALC.init(c(3,8,2))
result3<-alcove(parmSet,inp,exemplar,target3,n.train)
lines(result3\$accu,type='o',lwd=2,pch=20,col=3)
set.seed(seed.num)
parmSet<-ALC.init(c(3,8,2))
result4<-alcove(parmSet,inp,exemplar,target4,n.train)
lines(result4\$accu,type='o',lwd=2,pch=20,col=4)
set.seed(seed.num)
parmSet<-ALC.init(c(3,8,2))
result5<-alcove(parmSet,inp,exemplar,target5,n.train)
lines(result5\$accu,type='o',lwd=2,pch=20,col=5)
set.seed(seed.num)
parmSet<-ALC.init(c(3,8,2))
result6<-alcove(parmSet,inp,exemplar,target6,n.train)
lines(result6\$accu,type='o',lwd=2,pch=20,col=6)
legend("bottomright",paste("type",1:6,sep=''),col=1:6,lwd=2,pch=20)
### end ALCOVE

# evo. game
food = 6; cost = 10
A = 0.5*(food-cost); B = food
C = 0; D = food/2
pay.mat = matrix(c(A,B,C,D),nrow=2,byrow=TRUE)
dt = 0.01;max_time=1000
p = rep(0,max_time)
q = rep(0,max_time)
p[1] = 0.2
q[1] = 1 - p[1]
for (t in 1:max_time){
prob.mat = outer(c(p[t],q[t]),c(p[t],q[t]))
W.ave = sum(prob.mat*pay.mat)
W.h = sum(c(p[t],q[t])*pay.mat[1,])
W.d = sum(c(p[t],q[t])*pay.mat[2,])
p[(t+1)] = p[t]+(p[t]*(W.h-W.ave)/W.ave)*dt
q[(t+1)] = q[t]+(q[t]*(W.d-W.ave)/W.ave)*dt
}

plot(p,type='l',lwd=2,col='red',ylim=c(0,1))
lines(q,type='l',lwd=2,col='blue')

# cake game
n.cake = 10
pay.mat = matrix(0,n.cake+1,n.cake+1)
for (i.cake in 1:n.cake){
pay.mat[(i.cake+1),1:(n.cake-i.cake+1)] =i.cake
}

p.cake = runif(n.cake+1)
p.cake = p.cake/sum(p.cake)
max.time = 50
dt = 0.01
t = seq(0,max.time,dt)
n.iter = length(t)
p.hist = matrix(0,nrow = n.iter, ncol = (n.cake+1))
p.hist[1,] = p.cake

for (i.time in 2:n.iter){
W = colSums(p.cake*t(pay.mat))
W.ave = sum(outer(p.cake,p.cake)*pay.mat)
p.cake = p.cake + p.cake*(W - W.ave)/W.ave * dt
p.hist[i.time,] = p.cake
}

plot(p.hist[,1],ylim=c(0,1),type='l',lwd=2,ylab = 'Proportion',xlab="time")
for (i.strat in 2:(n.cake+1)){
lines(p.hist[,i.strat],col=i.strat,lwd=2)
}
legend("topleft",paste("request = ",0:10),col=1:(n.cake+1),lwd =2,cex=0.75)

### fighting couple
fighting_couple<-function(a,b,c,d) {
timeSep=0.05;ts=seq(1,50,timeSep);n_ts=length(ts)
x=matrix(0,nrow=n_ts,ncol=1)
y=matrix(0,nrow=n_ts,,ncol=1)
initX=c(rep(-40,5),rep(-20,5),rep(0,5),rep(20,5),rep(40,5))
initY=c(rep(c(-40,-20,0,20,40),5))
initX=initX[-13]
initY=initY[-13]
lengthINI=length(initX)
for (i_ini in 1:lengthINI) {
x[1]=initX[i_ini];y[1]=initY[i_ini];
for (i_gen in 2:n_ts) {
x[i_gen]=x[i_gen-1]+(a*x[i_gen-1]+b*y[i_gen-1])*timeSep
y[i_gen]=y[i_gen-1]+(c*x[i_gen-1]+d*y[i_gen-1])*timeSep
}
if (i_ini==1) {
plot(x,y,xlim=c(-50,50),ylim=c(-50,50),col=4,type='l',lwd=2,
xlab="X's Action",ylab="Y's Action")
arrows(x[2],y[2],x[3],y[3],col=4,lwd=2,length=0.15)} else {
lines(x, y, col=4, lwd=2)
arrows(x[2], y[2], x[3], y[3], col=4,lwd=2, length=0.15)
}
}
}

par(mfrow=c(2,3))
fighting_couple(-1,0.0,0.5,-1)
fighting_couple(-1,2,-1,-1)
fighting_couple(0,-1,1,0)
fighting_couple(1,-2,2,0)
fighting_couple(1,0,0.5,1)
fighting_couple(1,-4,-4,0)
```