クラスター分析入門　ファジィクラスタリングの理論と応用　宮本定明先生

```{r}

#FCM2 p.36

#エントロピー関数による正則化

data(iris)

n=150;c=3

X=iris[,!(colnames(iris) %in% "Species")]

X=matrix(unlist(X),nrow=nrow(X),ncol=ncol(X))*10

class=unique(iris$Species)

lam=0.1

ite=1000

uik=array(1,dim=c(n,c))

uik=matrix(sample(c(1:9),n*c,replace=T)/10,nrow=n)

for(l in 1:ite){

v=array(0,dim=c(length(class),ncol(X)))  
  
for(j in 1:c){  
    
x=X

vec=apply(x*c(uik[,j]),2,sum)/sum(uik[,j])

v[j,]=vec

}  
  
values=array(0,dim=dim(uik))

for(j in 1:c){
  
uik[,j]=exp(-lam*apply((t(t(X)-c(v[j,])))^2,1,sum)) 

values[,j]=apply((t(t(X)-c(v[j,])))^2,1,sum)

}

uik=uik/apply(uik,1,sum)

J=sum(uik*values)+sum(uik*log(uik))/lam

print(J)

}

```

```{r}


#FCM2 p.39

#2次正則化

#試行中

data(iris)

n=150;c=4

X=iris[,!(colnames(iris) %in% "Species")]

X=matrix(unlist(X),nrow=nrow(X),ncol=ncol(X))*10

class=unique(iris$Species)

lam=10^(-0)

ite=50

uik=array(1,dim=c(n,c))

uik=matrix(sample(c(1:9),n*c,replace=T)/10,nrow=n)

for(l in 1:ite){

v=array(0,dim=c(c,ncol(X)))  

for(j in 1:c){  

x=X

vec=apply(x*c(uik[,j]),2,sum)/sum(uik[,j])

v[j,]=vec

}  

values=array(0,dim=dim(uik))

for(j in 1:c){
  
values[,j]=apply((t(t(X)-c(v[j,])))^2,1,sum)

}

f=array(0,dim=dim(uik))

for(j in 1:c){

f[,j]=1+lam*(apply(values,1,sum)-c*values[,j])

}

uik=array(0,dim=dim(uik))

for(j in 1:c){
  
uik[,j]=(apply(ifelse(f>0,1,0)*values,1,sum)*lam+1)/apply(ifelse(f>0,1,0),1,sum)-lam*values[,j]  
  
}

uik=ifelse(uik>0,uik,0)

J=sum(uik*values)+sum(uik^2)/(2*lam)

print(J)

}



```


```{r}

#ファジィｃ多様体法 p.46

data(iris)

n=150;c=3

X=iris[,!(colnames(iris) %in% "Species")]

X=matrix(unlist(X),nrow=nrow(X),ncol=ncol(X))*10

class=unique(iris$Species)

Y=as.integer(iris$Species)*10

ite=100

uik=matrix(sample(c(1:9),n*c,replace=T)/10,nrow=n)

m=2;q=3

for(l in 1:ite){

w=array(0,dim=c(length(class),ncol(X)))  
  
for(j in 1:c){  
    
x=X

vec=apply(x*c(uik[,j]^m),2,sum)/sum(uik[,j]^m)

w[j,]=vec

}  


Dik=array(0,dim=c(c,n))

for(j in 1:c){

A=array(0,dim=c(ncol(X),ncol(X)))    

for(i in 1:nrow(X)){
  
A=A+(uik[i,j]^m)*t(t(c(X[i,]-w[j])))%*%t(c(X[i,]-w[j]))   
  
}  
  
s=eigen(A)$vectors[,1:q] 

Dik[j,]=apply(t(t(X)-c(w[j,]))^2,1,sum)-apply((t(t(X)-c(w[j,]))%*%s)^2,1,sum)

}

for(j in 1:c){
  
sub=(t(1/Dik)*c(Dik[j,]))^(1/(m-1))  
  
sub=apply(sub,1,sum)^(-1)

uik[,j]=sub

}  

J=sum(t(Dik)*uik^m)

print(J)

}

```

```{r}

#ファジィｃ回帰法 p.48

data(iris)

n=150;c=3

X=iris[,!(colnames(iris) %in% "Species")]

X=matrix(unlist(X),nrow=nrow(X),ncol=ncol(X))*10

class=unique(iris$Species)

Y=as.integer(iris$Species)*10

ite=100

uik=matrix(sample(c(1:9),n*c,replace=T)/10,nrow=n)

beta=array(1,dim=c(c,ncol(X)+1))

m=2

for(l in 1:ite){

Dik=array(0,dim=c(c,n))

for(j in 1:c){

Z=cbind(rep(1,nrow(X)),X)  
  
b=solve(t(Z)%*%diag(c(uik[,j]))%*%Z)%*%t(Z)%*%diag(c(uik[,j]))%*%Y  
  
beta[j,]=b

Dik[j,]=(Y-Z%*%b)^2
  
}

for(j in 1:c){
  
sub=(t(1/Dik)*c(Dik[j,]))^(1/(m-1))  
  
sub=apply(sub,1,sum)^(-1)

uik[,j]=sub

}  

J=sum(t(Dik)*uik^m)

print(J)

}

```


```{r}

#5.3 非類似度を用いたファジィクラスタリング　p.72

mat=matrix(0,nrow=10,ncol=10)

mat[1,]=c(0,8.89,10.16,15.84,16.84,33.01,32.69,43.47,37.78,39.31)
mat[2,]=c(0,0,9.64,15.37,20.45,30.75,34.09,41.61,39.26,41.47)
mat[3,]=c(0,0,0,6.12,11.2,22.98,24.56,33.61,29.75,31.86)
mat[4,]=c(0,0,0,0,8.48,17.20,18.75,27.64,23.91,26.23)
mat[5,]=c(0,0,0,0,0,20.77,16.25,29.52,21.17,22.49)
mat[6,]=c(0,0,0,0,0,0,13.01,10.86,15.61,19.06)
mat[7,]=c(0,0,0,0,0,0,0,15.99,5.19,7.69)
mat[8,]=c(rep(0,8),14.84,17.68)
mat[9,]=c(rep(0,9),3.45)

mat=mat+t(mat)

d=mat

c=3

uij=matrix(sample(1:9,c*nrow(d),replace=T),nrow=c)

w=0.1

eta=0.001

ite=10

for(l in 1:ite){
  
w_pre=w  
  
uij_pre=uij  
  
for(j in 1:c){
  
u_vec=uij[j,]  

for(i in 1:ncol(uij)){
  
uij[j,i]=uij[j,i]-eta*sum((u_vec[i]-u_vec)*w*(w*(u_vec[i]-u_vec)^2-d[i,]))  
  
}
  
}  

du1=array(0,dim=c(ncol(uij),ncol(uij)))
  
du2=array(0,dim=c(ncol(uij),ncol(uij)))


for(j in 1:c){
  
u_vec=uij[j,]  

du1=du1+d*(t(matrix(rep(u_vec,length(u_vec)),ncol=length(u_vec)))-u_vec)^2
  
du2=du2+(t(matrix(rep(u_vec,length(u_vec)),ncol=length(u_vec)))-u_vec)^4

}

w=sum(du1)/sum(du2)

cost=sum((uij-uij_pre)^2)+(w-w_pre)^2

print(cost)
  
}


```