hfc.m

The neuro-fuzzy software for identification and data analysis has been implemented in the MATLAB lan

function [centers,bases]=hfc(c,data)

%HFC		Hyperplane Fuzzy Clustering
%
%  [centers,bases]=hfc(c,data)
%
%	centers[no_rules,no_var]   centers of the membership functions
%	bases[no_rules,no_var]	   bases   of the membership functions   	
%
%	c			complexity (number of fuzzy rules)
%	data			the training set matrix
%
	
%___________________________________________________20 April 1996_______

global memb_fun


alpha=.8;
m=2;

[N,p]=size(data);

ranges=[min(data);max(data)];


u=zeros(c,N);
distance=zeros(c,N);
D=zeros(c,N);
gamma=zeros(c,N);

%Random initialization of u with crisp values
for k=1:N
	[maximum,pos]=max(rand(c,1));
	u(pos,k)=1;
end



first=1;
repeat=1;

while (repeat)
	
	for i=1:c
		
		S=zeros(p,p);			%Scatter Matrix
		hyperplane=zeros(p-1,p);	
						
		
		if sum(u(i,:).^m)
			v(i,:)=(u(i,:).^m)*data./sum(u(i,:).^m);
		else
			for j=1:p
				v(i,j)=ranges(1,j)+rand*(ranges(2,j)-ranges(1,j));
			end	
		end
		
		
		for k=1:N
			S=S+(u(i,k).^m)*(data(k,:)-v(i,:))'*(data(k,:)-v(i,:));
		end
	
		[V,Lambda]=eig(S);

		[lambda_min,pos_min]=min(diag(Lambda));

		V(:,pos_min)=[];
		hyperplane=V';
				
		if (trace(Lambda)-lambda_min)
			epsilon=1-lambda_min/((trace(Lambda)-lambda_min));
		else
			epsilon=0;
		end


		for k=1:N
			distance(i,k)=norm(data(k,:)-v(i,:));
			D(i,k)=distance(i,k)*sin(subspace((data(k,:)-v(i,:))',hyperplane'));
			gamma(i,k)=distance(i,k).^2-epsilon*D(i,k).^2;
		end
	end


	
	for i=1:c
		for k=1:N
			temp=0;
			for j=1:c
				if gamma(j,k)
					temp=temp+gamma(i,k)/gamma(j,k);
				end
			end
			
			if temp
				u(i,k)=temp^(-1/(m-1));
			else
				u(i,k)=0;
			end
		end
	end


	
	J=sum(sum((u.^m).*(alpha*D.^2+(1-alpha)*distance.^2)));
	
	if first
		first=0;
	else
		if (J>J_old-0.3)
			repeat=0;
		end
	end
	
	J_old=J;
end


centers=v(:,1:p-1);
bases=zeros(size(centers));
default_base=(ranges(2,1:p-1)-ranges(1,1:p-1))./20;


if strcmp(memb_fun,'gaussian')
	for i=1:c
		related=(u(i,:)>0.05);
		related_data=data(related,1:p-1);
		if size(related_data,1)>1
			bases(i,:)=(std(related_data)).^2;
		else
			bases(i,:)=default_base;
		end
	end
end

if strcmp(memb_fun,'triangular')
	for i=1:c
		related=(u(i,:)>0.05);
		related_data=data(related,1:p-1);
		no_related_data=size(related_data,1);
		
		bas=abs(ones(no_related_data,1)*centers(i,:)-related_data);
		
		if no_related_data>1
			bases(i,:)=2*max(bas);
		elseif no_related_data==1
			bases(i,:)=2*bas;
		else
			bases(i,:)=default_base;
		end
	end
end