% what is the difference from the previous Bayesian. -- it is a vote from different models in the previous files, 

% *** not sure whether I'm wasting my time to do sth which is not supported by theory. 
% There is a problem when you decompose them.
% don't worry, if it is the parameter learn from data -- then from the point of SLP, that's completely correct. 
% on the other hand, is it where the problems come from?
% the link matrix is decomposed into several rules with one-to-one mapping



testBayePrediction(EI,Accuracy):-
	ex(EI,Goal,PosSign), Goal=concentration(M,ActualChange,Time),
	bayesPredicting(Goal,HChange), %bayesPredicting(M,Time,HChange),
	write('----------------- Check predicted results'),write(EI),write(HChange),nl,
	(HChange=[ActualChange] ->
		Accuracy=1,
		write('Correctly predicted');
		(HChange=[PossChange1,PossChange2]->
			Accuracy=1/2,write('Partly right');
			Accuracy=0, write('NO :(')
		)		

		/*(member(ActualChange,HChange)->
			length(HChange,N),Accuracy=1/N,write('Partly right');
			Accuracy=0, write('NO :(')
		)*/
	).

	
/*	List=[up-1,down-1,no_change-1],
	max_of_3values(List,Keys-P),
	write(Keys-P).

*/


% this predicate can be turned into a more general predicate, which choose the two values (yes or no)
bayesPredicting(Goal,Changes):-
	Goal = concentration(M,ActualChange,Time),
	oneGoalSumP(concentration(M,up,Time),P_up),
	oneGoalSumP(concentration(M,down,Time),P_down),
	oneGoalSumP(concentration(M,no_change,Time),P_noChange),
	write([up-P_up, down-P_down, no_change-P_noChange]),
	max_of_3values([up-P_up, down-P_down, no_change-P_noChange],Changes-Probability).


oneGoalSumP(Goal,SumP):-
	testPoint(1),allSucceedPaths(Goal,Ts), testPoint(1),
	print_list(Ts),

	(Ts=[]->
		SumP=0;
		maplist(processP,Ts,TsP),
%processP(Ts,TsP), % in terms of computational efficiency, it is better to attached the probability -- or you can just record the probability--no care about the exact used clause. 
				% another problem -- you need to augment the hypothesized clauses with TSoFar,T in order to track proving paths. 
		sumList(TsP,SumP)
	).


processP([],1).
processP([Cla|T],T_Probability):-
	processP(T,PreT_Probability),
	p(Cla,Cla_Probability),
	T_Probability = Cla_Probability*PreT_Probability.
	




allSucceedPaths(concentration(MID,Change,Time),Ts):-
	set(depthLimit,DepthLimit),
	findall(OneTr,
		concentration(DepthLimit,[],OneTr,[],MID,Change,Time), %AugmentedGoal
		Ts).



% output the whole (both key and probability )
% a list of keys with the same max value
max_of_3values([],[]-0).
max_of_3values([Key-Probaility|List],Max):-
	max_of_3values(List,Keys-P),
	(Probaility > P ->
		Max=[Key]-Probaility;
		( Probaility = P ->
			Max=[Key|Keys]-P;
			Max=Keys-P 	% put it as before
		)
	).



claAveragePrediction(TsPs,FactsPs):-
	rb_new(Tree0),
	claAveragePredictionRecord(TsPs,Tree0,Tree),
	rb_visit(Tree,FactsPs0), write('OK after counting'),%write(TIsEIs),nl,
	maplist(normalizaStates,FactsPs0,FactsPs),
	print_list(FactsPs),nl.

claAveragePredictionRecord([],Tree,Tree).
claAveragePredictionRecord([P-T|All_T_Records],TreeSoFar,Tree):-
	eachModel(T,P,TreeSoFar,NTreeSoFar),
	claAveragePredictionRecord(All_T_Records,NTreeSoFar,Tree).

eachModel([],P,Tree,Tree).
eachModel([ClaI|ClaIs],P,TreeSoFar,Tree):-
	ClaI = [enzyme_state-{Score,EC,State,Time}],
	(rb_lookup(EC-Time,HChangeSoFar,TreeSoFar)->
		updataP(State,P,HChangeSoFar,NewHChange),
		rb_update(TreeSoFar,EC-Time,NewHChange,NTreeSoFar)
		;
		initialize(State,P,InitialHChange),	
		rb_insert(TreeSoFar,EC-Time,InitialHChange,NTreeSoFar)
	),	
	eachModel(ClaIs,P,NTreeSoFar,Tree).



initialize(activated,P,[activated-P,inhibited-0,no_change-0]).
initialize(inhibited,P,[activated-0,inhibited-P,no_change-0]).
initialize(no_change,P,[activated-0,inhibited-0,no_change-P]).
	
updataP(activated,P,[activated-PrePa,inhibited-PrePi,no_change-PrePno],[activated-NewPa,inhibited-NewPi,no_change-NewPno]):- NewPa is PrePa+P, NewPi=PrePi, NewPno=PrePno.
updataP(inhibited,P,[activated-PrePa,inhibited-PrePi,no_change-PrePno],[activated-NewPa,inhibited-NewPi,no_change-NewPno]):- NewPa=PrePa, NewPi is PrePi+P, NewPno=PrePno.
updataP(no_change,P,[activated-PrePa,inhibited-PrePi,no_change-PrePno],[activated-NewPa,inhibited-NewPi,no_change-NewPno]):- NewPa=PrePa, NewPi=PrePi, NewPno is PrePno+P.
	

normalizaStates(EC-Time-HChange,EC-Time-NewHChange):-
	maplist(getStateP,HChange,Ps),
	sumList(Ps,Total),
	
	HChange=[activated-PrePa,inhibited-PrePi,no_change-PrePno],
	NewPa is PrePa/Total,
	NewPi is PrePi/Total,
	NewPno is PrePno/Total,
	NewHChange=[activated-NewPa,inhibited-NewPi,no_change-NewPno].

getStateP(State-P,P).