コンジョイント分析

#数量化1類

data=data.frame(販売率=c(3.0,1.8,1.5,3.3,2.2,2.0,3.5,2.0,1.7,2.3),大都市=c(1,0,0,1,1,1,1,0,1,1),地方都市=c(0,1,1,0,0,0,0,1,0,0),温暖=c(0,1,0,0,0,0,1,1,0,0),普通=c(1,0,1,1,0,0,0,0,0,0),寒冷=c(0,0,0,0,1,1,0,0,1,1))

X=data[,colnames(data)!="温暖"];X=X[,2:ncol(X)]

x=array(0,dim=c(length(colnames(X)),length(colnames(X))))

for(i in 1:ncol(X)){
for(j in 1:ncol(X)){  


  x[i,j]=sum(X[,i]*X[,j])


}
}

xy=array(0,dim=c(ncol(X),1))

for(i in 1:ncol(X)){

xy[i]=sum(X[,i]*data$販売率) 

}

A=solve(x)%*%xy

y=c()

for(j in 1:nrow(X)){

y=c(y,sum(X[j,]*A))  

}

U=y



#MONANOVA

for(k in 2:length(y)){

if(U[k]<U[k-1]){

MU_data=data.frame(num=1:(k-1),MU=0,sig=0)

for(l in 1:(k-1)){  

MU=sum(U[k]+U[(k-1):(k-l)])/(l+1)

if(k-l-1==0){

MU_data$MU[l]=MU;MU_data$sig[l]=ifelse(MU>0,1,0)

}else{

MU_data$MU[l]=MU;MU_data$sig[l]=ifelse(MU>U[(k-l-1)],1,0)

}  

}

no=min(MU_data$num[MU_data$sig==1])

MU=MU_data$MU[no]

U[k:(k-no)]=MU

}

U_ave=c()  

for(l in 1:length(U)){

U_ave=c(U_ave,mean(U[1:l]))  

}  


eta=sqrt(sum((U-U_ave)^2)/sum((U-mean(U))^2))

print(eta)

}  

#効用推定法


X=as.matrix(data[,2:ncol(data)]);Y=c(data$販売率)

a=rep(1,ncol(X));b=rep(1,nrow(X))

ite=10000;delta=0.01

for(j in 1:ite){

a_pre=a;b_pre=b  

A=sum((U-U_ave)^2);B=sum((U-mean(U))^2)

eta=sqrt(sum((U-U_ave)^2)/sum((U-mean(U))^2))


for(i in 1:length(a)){

da=2*sum(X[,i]*(U-U_ave))

db=2*sum(X[,i]*(U-mean(U)))

de=(da*eta/A-db*eta/B)/2

a[i]=a[i]-delta*de

}  

for(i in 1:length(b)){

da=2*sum((U-U_ave)*Y)

db=2*sum(U-mean(U)*Y)

de=(da*eta/A-db*eta/B)/2

b[i]=b[i]-delta*de

} 

U=c(X%*%a)+c(b*Y)

U_ave=c()

for(l in 1:length(U)){

U_ave=c(U_ave,mean(U[1:l]))  

}

print(sum(abs(c(a_pre,b_pre)-c(a,b))))

}