112 lines
4.1 KiB
Plaintext
112 lines
4.1 KiB
Plaintext
|
%%%%
|
||
|
|
%%%% ABO blood type --- bloodABO.psm
|
||
|
|
%%%%
|
||
|
|
%%%% Copyright (C) 2004,2006,2008
|
||
|
|
%%%% Sato Laboratory, Dept. of Computer Science,
|
||
|
|
%%%% Tokyo Institute of Technology
|
||
|
|
|
||
|
|
%% ABO blood type consists of A, B, O and AB. They are observable
|
||
|
|
%% (phenotypes) and determined by a pair of blood type genes (geneotypes).
|
||
|
|
%% There are three ABO genes, namely a, b and o located on the 9th
|
||
|
|
%% chromosome of a human being. There are 6 geneotypes ({a,a},{a,b},{a,o},
|
||
|
|
%% {b,b},{b,o},{o,o}) and each determines a blood type. For example {a,b}
|
||
|
|
%% gives blood type AB etc. Our task is to estimate frequencies of ABO
|
||
|
|
%% genes from a random sample of ABO blood type, assuming random mate.
|
||
|
|
|
||
|
|
%%-------------------------------------
|
||
|
|
%% Quick start : sample session
|
||
|
|
%%
|
||
|
|
%% ?- prism(bloodABO),go,print_blood.
|
||
|
|
%% % Learn parameters from randomly generated
|
||
|
|
%% % 100 samples with A:B:O:AB = 38:22:31:9
|
||
|
|
%%
|
||
|
|
%% ?- sample(bloodtype(X)).
|
||
|
|
%% % Pick up a person with blood type X randomly
|
||
|
|
%% % acccording to the currrent parameter settings
|
||
|
|
%%
|
||
|
|
%% ?- get_samples(100,bloodtype(X),_Gs),countlist(_Gs,Cs).
|
||
|
|
%% % Pick up 100 persons and get the frequencies
|
||
|
|
%% % of their blood types
|
||
|
|
%%
|
||
|
|
%% ?- probf(bloodtype(ab),E),print_graph(E).
|
||
|
|
%% % Print all explanations for blooodtype(ab) in
|
||
|
|
%% % a compressed form
|
||
|
|
%%
|
||
|
|
%% ?- prob(bloodtype(ab),P).
|
||
|
|
%% % P is the probability of bloodtype(ab) being true
|
||
|
|
%%
|
||
|
|
%% ?- viterbif(bloodtype(ab)).
|
||
|
|
%% ?- viterbif(bloodtype(ab),P,E),print_graph(E).
|
||
|
|
%% ?- viterbi(bloodtype(ab),P).
|
||
|
|
%% % P is the probability of a most likely
|
||
|
|
%% % explanation E for bloodtype(ab).
|
||
|
|
%%
|
||
|
|
%% ?- viterbit(bloodtype(ab)).
|
||
|
|
%% % Print the most likely explanation for
|
||
|
|
%% % bloodtype(ab) in a tree form.
|
||
|
|
|
||
|
|
go:- learn_bloodtype(100).
|
||
|
|
|
||
|
|
%%-------------------------------------
|
||
|
|
%% Declarations:
|
||
|
|
|
||
|
|
:- set_prism_flag(data_source,file('bloodtype.dat')).
|
||
|
|
% When we run learn/0, the data are supplied
|
||
|
|
% by `bloodtype.dat'.
|
||
|
|
|
||
|
|
values(gene,[a,b,o],[0.5,0.2,0.3]).
|
||
|
|
% We declare msw(gene,V) s.t. V takes on
|
||
|
|
% one of the genes {a,b,o} when executed,
|
||
|
|
% with the freq.: a 50%, b 20%, o 30%.
|
||
|
|
|
||
|
|
%%------------------------------------
|
||
|
|
%% Modeling part:
|
||
|
|
|
||
|
|
bloodtype(P) :-
|
||
|
|
genotype(X,Y),
|
||
|
|
( X=Y -> P=X
|
||
|
|
; X=o -> P=Y
|
||
|
|
; Y=o -> P=X
|
||
|
|
; P=ab
|
||
|
|
).
|
||
|
|
|
||
|
|
genotype(X,Y) :- msw(gene,X),msw(gene,Y).
|
||
|
|
% We assume random mate. Note that msw(gene,X)
|
||
|
|
% and msw(gene,Y) are i.i.d. (independent and
|
||
|
|
% identically distributed) random variables
|
||
|
|
% in Prism because they have the same id but
|
||
|
|
% different subgoals.
|
||
|
|
|
||
|
|
%%------------------------------------
|
||
|
|
%% Utility part:
|
||
|
|
|
||
|
|
learn_bloodtype(N) :- % Learn parameters from N observations
|
||
|
|
random_set_seed(214857), % Set seed of the random number generator
|
||
|
|
gen_bloodtype(N,Gs),!, % Sample bloodtype/1 of size N
|
||
|
|
learn(Gs). % Perform search and graphical EM learning
|
||
|
|
% learn. % <= when using the file `bloodtype.dat'
|
||
|
|
|
||
|
|
gen_bloodtype(N,Gs) :-
|
||
|
|
N > 0,
|
||
|
|
random_select([a,b,o,ab],[0.38,0.22,0.31,0.09],X),
|
||
|
|
Gs = [bloodtype(X)|Gs1], % Sample a blood type with an empirical
|
||
|
|
N1 is N-1,!, % ratio for Japanese people.
|
||
|
|
gen_bloodtype(N1,Gs1).
|
||
|
|
gen_bloodtype(0,[]).
|
||
|
|
|
||
|
|
print_blood :-
|
||
|
|
prob(bloodtype(a),PA),prob(bloodtype(b),PB),
|
||
|
|
prob(bloodtype(o),PO),prob(bloodtype(ab),PAB),
|
||
|
|
nl,
|
||
|
|
format("P(A) = ~6f~n",[PA]),
|
||
|
|
format("P(B) = ~6f~n",[PB]),
|
||
|
|
format("P(O) = ~6f~n",[PO]),
|
||
|
|
format("P(AB) = ~6f~n",[PAB]).
|
||
|
|
|
||
|
|
print_gene :-
|
||
|
|
get_sw(gene,[_,[a,b,o],[GA,GB,GO]]),
|
||
|
|
nl,
|
||
|
|
format("P(a) = ~6f~n",[GA]),
|
||
|
|
format("P(b) = ~6f~n",[GB]),
|
||
|
|
format("P(o) = ~6f~n",[GO]).
|