function phenotypeFactor = phenotypeGivenGenotypeMendelianFactor(isDominant, genotypeVar, phenotypeVar)

 % This function computes the probability of each phenotype given the

 % different genotypes for a trait.  It assumes that there is 1 dominant

 % allele and 1 recessive allele.

 %

 % If you do not have much background in genetics, you should read the

 % on-line Appendix or watch the Khan Academy Introduction to Heredity Video

 % (http://www.khanacademy.org/video/introduction-to-heredity?playlist=Biology)

 % before you start coding.

 %

 % For the genotypes, assignment 1 maps to homozygous dominant, assignment 2

 % maps to heterozygous, and assignment 3 maps to homozygous recessive.  For

 % example, say that there is a gene with two alleles, dominant allele A and

 % recessive allele a.  Assignment 1 would map to AA, assignment 2 would

 % make to Aa, and assignment 3 would map to aa.  For the phenotypes,

 % assignment 1 maps to having the physical trait, and assignment 2 maps to

 % not having the physical trait.

 %

 % THE VARIABLE TO THE LEFT OF THE CONDITIONING BAR MUST BE THE FIRST

 % VARIABLE IN THE .var FIELD FOR GRADING PURPOSES

 %

 % Input:

 %   isDominant: 1 if the trait is caused by the dominant allele (trait

 %   would be caused by A in example above) and 0 if the trait is caused by

 %   the recessive allele (trait would be caused by a in the example above)

 %   genotypeVar: The variable number for the genotype variable (goes in the

 %   .var part of the factor)

 %   phenotypeVar: The variable number for the phenotype variable (goes in

 %   the .var part of the factor)

 %

 % Output:

 %   phenotypeFactor: Factor denoting the probability of having each

 %   phenotype for each genotype

 phenotypeFactor = struct('var', [], 'card', [], 'val', []);

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 %INSERT YOUR CODE HERE

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 % Fill in phenotypeFactor.var.  This should be a 1-D row vector.

 % Fill in phenotypeFactor.card.  This should be a 1-D row vector.

 phenotypeFactor.var = [phenotypeVar genotypeVar];

 phenotypeFactor.card= [2 3];

 array = IndexToAssignment(1:prod(phenotypeFactor.card),phenotypeFactor.card);

     if isDominant

         phenotypeFactor.val = ones(1, prod(phenotypeFactor.card));

         for i = 1:length(array)

             if array(i,2)==3&&array(i,1)==1

                 phenotypeFactor.val(i)=0;

             elseif ismember(array(i,2),[1 2])&&array(i,1)==2

                 phenotypeFactor.val(i)=0;

             end

         end

     else

         phenotypeFactor.val = zeros(1, prod(phenotypeFactor.card));

         for i = 1:length(array)

             if array(i,2)==3&&array(i,1)==1

                 phenotypeFactor.val(i)=1;

             end

         end

     end

 end

 % Replace the zeros in phentoypeFactor.val with the correct values.

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function phenotypeFactor = phenotypeGivenGenotypeFactor(alphaList, genotypeVar, phenotypeVar)

% This function computes the probability of each phenotype given the

% different genotypes for a trait. Genotypes (assignments to the genotype

% variable) are indexed from 1 to the number of genotypes, and the alphas

% are provided in the same order as the corresponding genotypes so that the

% alpha for genotype assignment i is alphaList(i).

%

% For the phenotypes, assignment 1 maps to having the physical trait, and

% assignment 2 maps to not having the physical trait.

%

% THE VARIABLE TO THE LEFT OF THE CONDITIONING BAR MUST BE THE FIRST

% VARIABLE IN THE .var FIELD FOR GRADING PURPOSES

%

% Input:

%   alphaList: Vector of alpha values for each genotype (n x 1 vector,

%   where n is the number of genotypes) -- the alpha value for a genotype

%   is the probability that a person with that genotype will have the

%   physical trait

%   genotypeVar: The variable number for the genotype variable (goes in the

%   .var part of the factor)

%   phenotypeVar: The variable number for the phenotype variable (goes in

%   the .var part of the factor)

%

% Output:

%   phenotypeFactor: Factor in which the val has the probability of having

%   each phenotype for each genotype combination (note that this is the

%   FULL CPD with no evidence observed)

phenotypeFactor = struct('var', [], 'card', [], 'val', []);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%INSERT YOUR CODE HERE

% The number of genotypes is the length of alphaList.

% The number of phenotypes is 2.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

% Fill in phenotypeFactor.var.  This should be a 1-D row vector.

% Fill in phenotypeFactor.card.  This should be a 1-D row vector.

phenotypeFactor.var = [phenotypeVar genotypeVar];

phenotypeFactor.card= [2 3];

phenotypeFactor.val = zeros(1, prod(phenotypeFactor.card));

% Replace the zeros in phentoypeFactor.val with the correct values.

assign_ = IndexToAssignment(1:prod(phenotypeFactor.card),[2 3]);

for i = 1:prod(phenotypeFactor.card)

    if assign_(i,1) == 1

        switch assign_(i,2)

            case 1

                phenotypeFactor.val(i) = alphaList(1);

            case 2

                phenotypeFactor.val(i) = alphaList(2);

            case 3

                phenotypeFactor.val(i) = alphaList(3);

            otherwise

                'mistake in genotype'

                return

        end

    else

         switch assign_(i,2)

            case 1

                phenotypeFactor.val(i) = 1-alphaList(1);

            case 2

                phenotypeFactor.val(i) = 1-alphaList(2);

            case 3

                phenotypeFactor.val(i) = 1-alphaList(3);

            otherwise

                'mistake in genotype'

                return

        end

    end

end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 function genotypeFactor = genotypeGivenAlleleFreqsFactor(alleleFreqs, genotypeVar)

 % alleleFreqs = [0.1; 0.9];

 % genotypeVar = 1;

 % This function computes the probability of each genotype given the allele

 % frequencies in the population.

 % Note that we assume that the copies of the gene are independent.  Thus,

 % knowing the allele for one copy of the gene does not affect the

 % probability of having each allele for the other copy of the gene.  As a

 % result, the probability of a genotype is the product of the frequencies

 % of its constituent alleles (or twice that product for heterozygous

 % genotypes).

 % Input:

 %   alleleFreqs: An n x 1 vector of the frequencies of the alleles in the

 %   population, where n is the number of alleles

 %   genotypeVar: The variable number for the genotype (goes in the .var

 %   part of the factor)

 %

 % Output:

 %   genotypeFactor: Factor in which the val has the probability of having

 %   each genotype (note that this is the FULL CPD with no evidence

 %   observed)

 % The number of genotypes is (number of alleles choose 2) + number of

 % alleles -- need to add number of alleles at the end to account for

 % homozygotes

 genotypeFactor = struct('var', [], 'card', [], 'val', []);

 numAlleles = length(alleleFreqs);

 % Each allele has an ID that is the index of its allele frequency in the

 % allele frequency list.  Each genotype also has an ID.  We need allele and

 % genotype IDs so that we know what genotype and alleles correspond to each

 % probability in the .val part of the factor.  For example, the first entry

 % in .val corresponds to the probability of having the genotype with

 % genotype ID 1, which consists of having two copies of the allele with

 % allele ID 1.  There is a mapping from a pair of allele IDs to genotype

 % IDs and from genotype IDs to a pair of allele IDs below; we compute this

 % mapping using generateAlleleGenotypeMappers(numAlleles). (A genotype

 % consists of 2 alleles.)

 [allelesToGenotypes, genotypesToAlleles] = generateAlleleGenotypeMappers(numAlleles);

 % One or both of these matrices might be useful.

 %

 %   1.  allelesToGenotypes: n x n matrix that maps pairs of allele IDs to

 %   genotype IDs, where n is the number of alleles -- if

 %   allelesToGenotypes(i, j) = k, then the genotype with ID k comprises of

 %   the alleles with IDs i and j

 %

 %   2.  genotypesToAlleles: m x 2 matrix of allele IDs, where m is the

 %   number of genotypes -- if genotypesToAlleles(k, :) = [i, j], then the

 %   genotype with ID k is comprised of the allele with ID i and the allele

 %   with ID j

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 %INSERT YOUR CODE HERE

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 % Fill in genotypeFactor.var.  This should be a 1-D row vector.

 % Fill in genotypeFactor.card.  This should be a 1-D row vector.

 genotypeFactor.var = genotypeVar;

 genotypeFactor.card = length(genotypesToAlleles(:,1));

 genotypeFactor.val = zeros(1, prod(genotypeFactor.card));

 % Replace the zeros in genotypeFactor.val with the correct values.

 n_ = length(alleleFreqs);

 for i = 1:n_

     for j = 1:n_

         num_of_val = allelesToGenotypes(i,j);

         genotypeFactor.val(num_of_val) = genotypeFactor.val(num_of_val)+alleleFreqs(i)*alleleFreqs(j) ;

     end

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 function genotypeFactor = genotypeGivenParentsGenotypesFactor(numAlleles, genotypeVarChild, genotypeVarParentOne, genotypeVarParentTwo)

 % This function computes a factor representing the CPD for the genotype of

 % a child given the parents' genotypes.

 % THE VARIABLE TO THE LEFT OF THE CONDITIONING BAR MUST BE THE FIRST

 % VARIABLE IN THE .var FIELD FOR GRADING PURPOSES

 % When writing this function, make sure to consider all possible genotypes

 % from both parents and all possible genotypes for the child.

 % Input:

 %   numAlleles: int that is the number of alleles

 %   genotypeVarChild: Variable number corresponding to the variable for the

 %   child's genotype (goes in the .var part of the factor)

 %   genotypeVarParentOne: Variable number corresponding to the variable for

 %   the first parent's genotype (goes in the .var part of the factor)

 %   genotypeVarParentTwo: Variable number corresponding to the variable for

 %   the second parent's genotype (goes in the .var part of the factor)

 %

 % Output:

 %   genotypeFactor: Factor in which val is probability of the child having

 %   each genotype (note that this is the FULL CPD with no evidence

 %   observed)

 % The number of genotypes is (number of alleles choose 2) + number of

 % alleles -- need to add number of alleles at the end to account for homozygotes

 genotypeFactor = struct('var', [], 'card', [], 'val', []);

 % Each allele has an ID.  Each genotype also has an ID.  We need allele and

 % genotype IDs so that we know what genotype and alleles correspond to each

 % probability in the .val part of the factor.  For example, the first entry

 % in .val corresponds to the probability of having the genotype with

 % genotype ID 1, which consists of having two copies of the allele with

 % allele ID 1, given that both parents also have the genotype with genotype

 % ID 1.  There is a mapping from a pair of allele IDs to genotype IDs and

 % from genotype IDs to a pair of allele IDs below; we compute this mapping

 % using generateAlleleGenotypeMappers(numAlleles). (A genotype consists of

 % 2 alleles.)

 [allelesToGenotypes, genotypesToAlleles] = generateAlleleGenotypeMappers(numAlleles);

 % One or both of these matrices might be useful.

 %

 %   1.  allelesToGenotypes: n x n matrix that maps pairs of allele IDs to

 %   genotype IDs, where n is the number of alleles -- if

 %   allelesToGenotypes(i, j) = k, then the genotype with ID k comprises of

 %   the alleles with IDs i and j

 %

 %   2.  genotypesToAlleles: m x 2 matrix of allele IDs, where m is the

 %   number of genotypes -- if genotypesToAlleles(k, :) = [i, j], then the

 %   genotype with ID k is comprised of the allele with ID i and the allele

 %   with ID j

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 %INSERT YOUR CODE HERE

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 % Fill in genotypeFactor.var.  This should be a 1-D row vector.

 % Fill in genotypeFactor.card.  This should be a 1-D row vector.

 genotypeFactor.var = [genotypeVarChild, genotypeVarParentOne, genotypeVarParentTwo];

 var_num_ = length(genotypeFactor.var);

 var_length_ = length(genotypesToAlleles(:,1));

 genotypeFactor.card = linspace(var_length_,var_length_,var_num_);

 genotypeFactor.val = zeros(1, prod(genotypeFactor.card));

 assign_ = IndexToAssignment(1:prod(genotypeFactor.card),genotypeFactor.card);

 % Replace the zeros in genotypeFactor.val with the correct values.

 for i = 1:length(assign_)

     GeType_child = assign_(i,1);

     GeType_par1  = assign_(i,2);

     GeType_par2  = assign_(i,3);

     AllType_child = genotypesToAlleles(GeType_child,:);

     AllType_par1 = genotypesToAlleles(GeType_par1,:);

     AllType_par2 = genotypesToAlleles(GeType_par2,:);

     combination = [];

     length_tmp_allpar = length(AllType_par1);

     for j = 1:length_tmp_allpar

         for k = 1:length_tmp_allpar

             combination = [combination;AllType_par1(j) AllType_par2(k)];

         end

     end

     num_count = 0;

     length_tmp_combination = length(combination);

     for j = 1:length_tmp_combination

         combin_row = combination(j,:);

         if sort(AllType_child) == sort(combin_row)

             num_count = num_count+1;

         end

     end

     genotypeFactor.val(i)= num_count/length_tmp_combination;

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

  function phenotypeFactor = phenotypeGivenCopiesFactor(alphaList, numAlleles, geneCopyVarOne, geneCopyVarTwo, phenotypeVar)

 % alphaList = [0.8; 0.6; 0.1; 0.5; 0.05; 0.01];

 % numAlleles = 3;

 % geneCopyVarOne= 1;

 % geneCopyVarTwo = 2;

 % phenotypeVar = 3;

 % This function makes a factor whose values are the probabilities of

 % a phenotype given an allele combination.  Note that a person has one

 % copy of the gene from each parent.

 % In the factor, each assignment maps to the allele at the corresponding

 % location on the allele list, so allele assignment 1 maps to

 % alleleList{1}, allele assignment 2 maps to alleleList{2}, ....  For the

 % phenotypes, assignment 1 maps to having the physical trait, and

 % assignment 2 maps to not having the physical trait.

 % THE VARIABLE TO THE LEFT OF THE CONDITIONING BAR MUST BE THE FIRST

 % VARIABLE IN THE .var FIELD FOR GRADING PURPOSES

 % Input:

 %   alphaList: m x 1 vector of alpha values for the different genotypes,

 %   where m is the number of genotypes -- the alpha value for a genotype

 %   is the probability that a person with that genotype will have the

 %   physical trait

 %   numAlleles: int that is the number of alleles

 %   geneCopyVarOne: Variable number corresponding to the variable for

 %   the first copy of the gene (goes in the .var part of the factor)

 %   geneCopyVarTwo: Variable number corresponding to the variable for

 %   the second copy of the gene (goes in the .var part of the factor)

 %   phenotypeVar: Variable number corresponding to the variable for the

 %   phenotype (goes in the .var part of the factor)

 %

 % Output:

 %   phenotypeFactor: Factor in which the values are the probabilities of

 %   having each phenotype for each allele combination (note that this is

 %   the FULL CPD with no evidence observed)

 phenotypeFactor = struct('var', [], 'card', [], 'val', []);

 % Each allele has an ID. Each genotype also has an ID, which is the index

 % of its alpha in the list of alphas.  There is a mapping from a pair of

 % allele IDs to genotype IDs and from genotype IDs to a pair of allele IDs

 % below; we compute this mapping using

 % generateAlleleGenotypeMappers(numAlleles). (A genotype consists of 2

 % alleles.)

 [allelesToGenotypes, genotypesToAlleles] = generateAlleleGenotypeMappers(numAlleles);

 % One or both of these matrices might be useful.

 %

 %   1.  allelesToGenotypes: n x n matrix that maps pairs of allele IDs to

 %   genotype IDs, where n is the number of alleles -- if

 %   allelesToGenotypes(i, j) = k, then the genotype with ID k comprises of

 %   the alleles with IDs i and j

 %

 %   2.  genotypesToAlleles: m x 2 matrix of allele IDs, where m is the

 %   number of genotypes -- if genotypesToAlleles(k, :) = [i, j], then the

 %   genotype with ID k is comprised of the allele with ID i and the allele

 %   with ID j

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 %INSERT YOUR CODE HERE

 % The number of phenotypes is 2

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 % Fill in phenotypeFactor.var.  This should be a 1-D row vector.

 % Fill in phenotypeFactor.card.  This should be a 1-D row vector.

 phenotypeFactor.var = [phenotypeVar geneCopyVarOne geneCopyVarTwo];

 phenotypeFactor.card= [2 numAlleles numAlleles];

 phenotypeFactor.val = zeros(1, prod(phenotypeFactor.card));

 % Replace the zeros in phentoypeFactor.val with the correct values.

 assign_alle = IndexToAssignment(1:prod(phenotypeFactor.card),phenotypeFactor.card);

 for i = 1:prod(phenotypeFactor.card)

     alle_type = assign_alle(i,2:3);

         if assign_alle(i,1) == 1

             phenotypeFactor.val(i) = alphaList(allelesToGenotypes(alle_type(1),alle_type(2)));

         else

             phenotypeFactor.val(i) = 1 - alphaList(allelesToGenotypes(alle_type(1),alle_type(2)));

         end

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

  function factorList = constructGeneticNetwork(pedigree, alleleFreqs, alphaList)

 % pedigree = struct('parents', [0,0;1,3;0,0]);

 % pedigree.names = {'Ira','James','Robin'};

 % alleleFreqs = [0.1; 0.9];

 % alphaList = [0.8; 0.6; 0.1];

 % This function constructs a Bayesian network for genetic inheritance.  It

 % assumes that there are only 2 phenotypes.  It also assumes that either

 % both parents are specified or neither parent is specified.

 %

 % In Matlab, each variable will have a number.  We need a consistent way of

 % numbering variables when instantiating CPDs so that we know what

 % variables are involved in each CPD.  For example, if IraGenotype is in

 % multiple CPDs and IraGenotype is number variable 1 in a CPD, then

 % IraGenotype should be numbered as variable 1 in all CPDs and no other

 % variables should be numbered 1; thus, every time variable 1 appears in

 % our network, we will know that it refers to Ira's genotype.

 %

 % Here is how the variables should be numbered, for a pedigree with n

 % people:

 %

 % 1.  The first n variables should be the genotype variables and should

 % be numbered according to the index of the corresponding person in

 % pedigree.names; the ith person with name pedigree.names{i} has genotype

 % variable number i.

 %

 % 2.  The next n variables should be the phenotype variables and should be

 % numbered according to the index of the corresponding person in

 % pedigree.names; the ith person with name pedigree.names{i} has phenotype

 % variable number n+i.

 %

 % Here is an example of how the variable numbering should work: if

 % pedigree.names = {'Ira', 'James', 'Robin'} and pedigree.parents = [0, 0;

 % 1, 3; 0, 0], then the variable numbering is as follows:

 %

 % Variable 1: IraGenotype

 % Variable 2: JamesGenotype

 % Variable 3: RobinGenotype

 % Variable 4: IraPhenotype

 % Variable 5: JamesPhenotype

 % Variable 6: RobinPhenotype

 %

 % Input:

 %   pedigree: Data structure that includes names, genders, and parent-child

 %   relationships

 %   alleleFreqs: Frequencies of alleles in the population

 %   alphaList: m x 1 vector of alpha values for genotypes, where m is the

 %   number of genotypes -- the alpha value for a genotype is the

 %   probability that a person with that genotype will have the

 %   physical trait

 %

 % Output:

 %   factorList: Struct array of factors for the genetic network (In each

 %   factor, .var, .card, and .val are all row 1-D vectors.)

 numPeople = length(pedigree.names);

 % Initialize factors

 % The number of factors is twice the number of people because there is a

 % factor for each person's genotype and a separate factor for each person's

 % phenotype.  Note that the order of the factors in the list does not

 % matter.

 factorList(2*numPeople) = struct('var', [], 'card', [], 'val', []);

 numAlleles = length(alleleFreqs); % Number of alleles

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 %INSERT YOUR CODE HERE

 % Variable numbers:

 % 1 - numPeople: genotype variables

 % numPeople+1 - 2*numPeople: phenotype variables

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 

 for i = 1 : 2*numPeople

     if i <= numPeople

         pedigree_parents = pedigree.parents(i,:);

         if max(pedigree_parents) == 0

             tmp_ = genotypeGivenAlleleFreqsFactor(alleleFreqs,i);

         else

             tmp_ = genotypeGivenParentsGenotypesFactor(length(alleleFreqs),i,pedigree_parents(1),pedigree_parents(2));

         end

     else

         tmp_ = phenotypeGivenGenotypeFactor(alphaList,i-numPeople,i);

     end

     factorList(i) = tmp_;

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 function factorList = constructDecoupledGeneticNetwork(pedigree, alleleFreqs, alphaList)

 % pedigree = struct('parents', [0,0;1,3;0,0]);

 % pedigree.names = {'Ira','James','Robin'};

 % alleleFreqs = [0.1; 0.7; 0.2];

 % alleleListThree = {'F', 'f', 'n'};

 % alphaList = [0.8; 0.6; 0.1; 0.5; 0.05; 0.01];

 % This function constructs a Bayesian network for genetic inheritance.  It

 % assumes that there are only 2 phenotypes.  It also assumes that, in the

 % pedigree, either both parents are specified or neither parent is

 % specified.

 %

 % In Matlab, each variable will have a number.  We need a consistent way of

 % numbering variables when instantiating CPDs so that we know what

 % variables are involved in each CPD.  For example, if IraGenotype is in

 % multiple CPDs and IraGenotype is number variable 1 in a CPD, then

 % IraGenotype should be numbered as variable 1 in all CPDs and no other

 % variables should be numbered 1; thus, every time variable 1 appears in

 % our network, we will know that it refers to Ira's genotype.

 %

 % Here is how the variables should be numbered, for a pedigree with n

 % people:

 %

 % 1.  The first n variables should be the gene copy 1 variables and should

 % be numbered according to the index of the corresponding person in

 % pedigree.names; the ith person with name pedigree.names{i} has gene copy

 % 1 variable number i.  If the parents are specified, then gene copy 1 is the

 % copy inherited from pedigree.names{pedigree.parents(i, 1)}.

 %

 % 2.  The next n variables should be the gene copy 2 variables and should

 % be numbered according to the index of the corresponding person in

 % pedigree.names; the ith person with name pedigree.names{i} has gene copy

 % 2 variable number n+i.  If the parents are specified, then gene copy 2 is the

 % copy inherited from pedigree.parents(i, 2).

 %

 % 3.  The final n variables should be the phenotype variables and should be

 % numbered according to the index of the corresponding person in

 % pedigree.names; the ith person with name pedigree.names{i} has phenotype

 % variable number 2n+i.

 %

 % Here is an example of how the variable numbering should work: if

 % pedigree.names = {'Ira', 'James', 'Robin'} and pedigree.parents = [0, 0;

 % 1, 3; 0, 0], then the variable numbering is as follows:

 %

 % Variable 1: IraGeneCopy1

 % Variable 2: JamesGeneCopy1

 % Variable 3: RobinGeneCopy1

 % Variable 4: IraGeneCopy2

 % Variable 5: JamesGeneCopy2

 % Variable 6: RobinGeneCopy2

 % Variable 7: IraPhenotype

 % Variable 8: JamesPhenotype

 % Variable 9: RobinPhenotype

 %

 % Input:

 %   pedigree: Data structure that includes the names and parents of each

 %   person

 %   alleleFreqs: n x 1 vector of allele frequencies in the population,

 %   where n is the number of alleles

 %   alphaList: m x 1 vector of alphas for different genotypes, where m is

 %   the number of genotypes -- the alpha value for a genotype is the

 %   probability that a person with that genotype will have the physical

 %   trait

 %

 % Output:

 %   factorList: struct array of factors for the genetic network (In each

 %   factor, .var, .card, and .val are row 1-D vectors.)

 numPeople = length(pedigree.names);

 % Initialize factors

 % We Need 3*numPeople factors because, for each person, there is a factor

 % for the CPD at each of 2 parental copies of the gene and a factor for the

 % CPD at the phenotype Note that the order of the factors in the list does

 % not matter.

 factorList(3*numPeople) = struct('var', [], 'card', [], 'val', []);

 % Initialize factors

 factorList(3*numPeople) = struct('var', [], 'card', [], 'val', []);

 numAlleles = length(alleleFreqs); % Number of alleles

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 % INSERT YOUR CODE HERE

 % Variable numbers:

 % 1 - numPeople: first parent copy of gene variables

 % numPeople+1 - 2*numPeople: second parent copy of gene variables

 % 2*numPeople+1 - 3*numPeople: phenotype variables

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 for k = 1:3

     switch k

         case 1

             for i = 1:numPeople

                 pedigree_parents = pedigree.parents(i,:);

                 if max(pedigree_parents)==0

                     factorList(i) = childCopyGivenFreqsFactor(alleleFreqs,i);

                 else

                     factorList(i) = childCopyGivenParentalsFactor(numAlleles,i,pedigree_parents(1),pedigree_parents(2));

                 end

             end

         case 2

             for i = (numPeople+1):2*numPeople

                 pedigree_parents = pedigree.parents(i-numPeople,:);

                 if max(pedigree_parents)==0

                     factorList(i) = childCopyGivenFreqsFactor(alleleFreqs,i);

                 else

                     factorList(i) = childCopyGivenParentalsFactor(numAlleles,i,pedigree_parents(1),pedigree_parents(2));

                 end

             end

         case 3

             for i = 2*numPeople+1:3*numPeople

                  factorList(i) = phenotypeGivenCopiesFactor(alphaList,numAlleles,i,i-numPeople,i-2*numPeople);

             end

     end

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

  function phenotypeFactor = constructSigmoidPhenotypeFactor(alleleWeights, geneCopyVarOneList, geneCopyVarTwoList, phenotypeVar)

 % alleleWeights = {[3, -3], [0.9, -0.8]};

 % geneCopyVarOneList = [1; 2];

 % geneCopyVarTwoList = [4; 5];

 % phenotypeVar = 3;

 % This function takes a cell array of alleles' weights and constructs a

 % factor expressing a sigmoid CPD.

 %

 % You can assume that there are only 2 genes involved in the CPD.

 %

 % In the factor, for each gene, each allele assignment maps to the allele

 % whose weight is at the corresponding location.  For example, for gene 1,

 % allele assignment 1 maps to the allele whose weight is at

 % alleleWeights{1}(1) (same as w_1^1), allele assignment 2 maps to the

 % allele whose weight is at alleleWeights{1}(2) (same as w_2^1),....

 %

 % You may assume that there are 2 possible phenotypes.

 % For the phenotypes, assignment 1 maps to having the physical trait, and

 % assignment 2 maps to not having the physical trait.

 %

 % THE VARIABLE TO THE LEFT OF THE CONDITIONING BAR MUST BE THE FIRST

 % VARIABLE IN THE .var FIELD FOR GRADING PURPOSES

 %

 % Input:

 %   alleleWeights: Cell array of weights, where each entry is an 1 x n

 %   of weights for the alleles for a gene (n is the number of alleles for

 %   the gene)

 %   geneCopyVarOneList: m x 1 vector (m is the number of genes) of variable

 %   numbers that are the variable numbers for each of the first parent's

 %   copy of each gene (numbers in this list go in the .var part of the

 %   factor)

 %   geneCopyVarTwoList: m x 1 vector (m is the number of genes) of variable

 %   numbers that are the variable numbers for each of the second parent's

 %   copy of each gene (numbers in this list go in the .var part of the

 %   factor) -- Note that both copies of each gene are from the same person,

 %   but each copy originally came from a different parent

 %   phenotypeVar: Variable number corresponding to the variable for the

 %   phenotype (goes in the .var part of the factor)

 %

 % Output:

 %   phenotypeFactor: Factor in which the values are the probabilities of

 %   having each phenotype for each allele combination (note that this is

 %   the FULL CPD with no evidence observed)

 phenotypeFactor = struct('var', [], 'card', [], 'val', []);

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

 % INSERT YOUR CODE HERE

 % Note that computeSigmoid.m will be useful for this function.

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  

 % Fill in phenotypeFactor.var.  This should be a 1-D row vector.

 % Fill in phenotypeFactor.card.  This should be a 1-D row vector.

 phenotypeFactor.var = [phenotypeVar geneCopyVarOneList' geneCopyVarTwoList'];

 %num_of_alle_in_each_geno

 geno_var_card = [];

 for i=1:length(geneCopyVarOneList)

     geno_var_card = [geno_var_card max(size(alleleWeights{i}))];

 end

 phenotypeFactor.card =[2 geno_var_card geno_var_card];

 phenotypeFactor.val = zeros(1, prod(phenotypeFactor.card));

 % Replace the zeros in phentoypeFactor.val with the correct values.

 assign_allele = IndexToAssignment(1:prod(phenotypeFactor.card),phenotypeFactor.card);

 for i = 1:prod(phenotypeFactor.card)

     pheno_allele_ = assign_allele(i,:);

 % start from the second var

     k = 2;

     weight = 0;

 % a pair of gene come from *two* people

     for j = 1:2

 % get num of allele from one people

         for l = 1:length(geneCopyVarOneList)

             weight = weight+alleleWeights{l}(pheno_allele_(k));

             k = k+1;

         end

     end

     switch pheno_allele_(1)

         case 1

             phenotypeFactor.val(i) = computeSigmoid(weight);

         case 2

             phenotypeFactor.val(i) = 1-computeSigmoid(weight);

     end

 end

 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%