;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
;                                                                     ;
;                            cooper.lisp                              ;
;                                                                     ;
;                         by Darren Boulton                           ;
;                                                                     ;
;                      last updated: 1/10/2003                        ;
;                                                                     ;
;        Code to calculate the Cooper and Herskovists formula         ;
;                                                                     ;
;                      call (cooper model data prior)                       ;
;                                                                     ;
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;

;A routine to simplify reloading of the code when modifying 
(defun l ()
  (load "cooper.lisp"))

;A sample model used in the Cooper and Herskovits paper.
;This model is referred to as Bs1 in the paper.
;A model comprises a list of lists, each inner list containing
;the index of each parent of that node.
;ie ((parents of node 0) (parents of node 1) ... (parents of node N))
;In the model below:
;    node 0 has no parents
;    node 1 has only node 0 as a parent
;    node 2 has only node 1 as a parent
;This gives the chain: 0->1->2
(defun model1 ()
  '(() (0) (1)))

;A sample model used in the Cooper and Herskovits paper.
;This model is referred to as Bs2 in the paper.
;A model comprises a list of lists, each inner list containing
;the index of each parent of that node.
;ie ((parents of node 0) (parents of node 1) ... (parents of node N))
;In the model below:
;    node 0 has no parents
;    node 1 has only node 0 as a parent
;    node 2 has only node 0 as a parent
;This gives the v-structure:   1<-0->2
(defun model2 ()
  '(() (0) (0)))

;A sample data set used in the Cooper and Herskovits paper.
;A data set is a list of lists, each inner list is a sample
;of the joint space with no missing values.
(defun data ()
  '((present absent absent)
    (present present present)
    (absent absent present)
    (present present present)
    (absent absent absent)
    (absent present present)
    (present present present)
    (absent absent absent)
    (present present present)
    (absent absent absent)))

;Returns the factorial of i 
(defun factorial (i)
  ;if i < 2 then simply return 1
  (if (< i 2)
      1
    ;else return i multiplied by the factorial of i-1
    (* i (factorial (- i 1)))))

;Returns the number of nodes implied by the data
;This is equal to the number of items in a sample
;Example: (nodes (data)) returns 3
(defun nodes (data)
  (length (car data)))

;Returns a list containing the indices of the parents of node
;ie (parent1 parent2 ... parentN)
;Example: (parents 2 (model1)) returns (1)
(defun parents (node model)
  (nth node model))

;Returns a list of all states for a particular node
;Note that duplicates are included
;ie (states 0 (data)) returns (present present absent present ...
;If duplicates are to be removed us (unique (states 0 (data)))
;which will return (present absent)
;See the description of unique below.
(defun states (node data)
  (if (null data)
      nil
    (cons (nth node (car data)) (states node (cdr data)))))

;Returns the list with duplicates removed
;Note that the order of the list in not guarenteed
;ie (unique '(a b a)) returns (a b)
(defun unique (l)
  (if (null l)
      nil
    ;if the first element in the list is also present
    ;in the remainder of the list...
    (if (find (car l) (cdr l) :test #'equal)
	;then ignore it
	(unique (cdr l))
      ;else keep it at the start of the list
      (cons (car l) (unique (cdr l))))))

;Returns the value Nij as used in the Cooper and Herskovits paper
;Also represented as Skj in the BAI textbook
;Calculates the number of samples in the data set in which a 
;particular parent instantiation occurs.
;    parent_instantiation is the specified value of a set of parents (k)
(defun N (parent_instantiation data)
  (if (null data)
      0
    ;if this sample contains the parent instantiation...
    (if (parents_match parent_instantiation (car data))
	;then add one to the tally
	(+ 1 (N parent_instantiation (cdr data)))
      ;else don't add one to the tally
      (N parent_instantiation (cdr data)))))

;Returns the value Alpha-kjl in the C&H formula
;   node is the index of the current node (k)
;   state is the specified state (l)
;   parent_instantiation is the specified values of the node's parents (j)
(defun alpha (node state parent_instantiation data)
  (if (null data)
      0
    ;if the sample has the specified value for the node...
    (if (null (and (equal (nth node (car data)) state)
		   ;and contains the correct parent instantiation...
		   (parents_match parent_instantiation (car data))))
	;then add one to the tally
	(alpha node state parent_instantiation (cdr data))
      ;else don't add one to the tally
      (+ 1 (alpha node state parent_instantiation (cdr data))))))

;Returns true if a sample contains the specified parent instantiation
(defun parents_match (parent_instantiation sample)
  (if (null parent_instantiation)
      t
    (if (or (null (car parent_instantiation)) (equal (car parent_instantiation) (car sample)))
	(parents_match (cdr parent_instantiation) (cdr sample))
      nil)))

;Returns a list of parent instantiations in a specified data set for a specified
;node in a specified model.
;This routine replaces the values for all variables that are not parents
;of the node with nil. This allows samples to be matched for a particular parent
;instantiation by only comparing the non-nil entries.
;Example: (parents 2 (model2)) returns (0)
;         Therefore: (parent_instantiations 2 (model2) (data)) returns
;                                       ((absent nil nil) (present nil nil))
;A special case occurs when a node has no parents, which would return ((nil nil nil))
;This case is handled individually to just return nil, representing the empty set
(defun parent_instantiations (node model data)
  ;if there are no parent instantiations...
  (if (all_nil (car (unique (filter_data (parents node model) data))))
      ;then just return nil
      nil
    ;otherwise return the list
    (unique (filter_data (parents node model) data))))


;Returns t only if all entries in the list are nil
;This routine is only used to detect the case when there are
;no parent instantiations in a list (see above)
(defun all_nil (l)
  (if (null l)
      t
    (if (null (car l))
	(all_nil (cdr l))
      nil)))

;Replaces the values of all variables not included in the list of nodes with nil
;Operates on whole data sets
;Example: (filter_data '(0 2) '((absent absent present))) returns (absent nil present)
(defun filter_data (nodes data)
  (if (null data)
      nil
    ;simply call the filter_sample routine for each sample in the data set
    (cons (filter_sample 0 nodes (car data)) (filter_data nodes (cdr data)))))

;Replaces the values of all variables not included in the list of nodes with nil
;Operates on individual samples
;Example: (filter_data '(0 2) '((absent absent present))) returns (absent nil present)
;count is always called with the value 0 and is used to iterate through the nodes
(defun filter_sample (count nodes sample)
  (if (null sample)
      nil
    ;if this node is present in the list of nodes (we don't want it replaced)...
    (if (find count nodes :test #'equal)
	;then keep the original value
	(cons (car sample) (filter_sample (+ 1 count) nodes (cdr sample)))
      ;else replace the original value with nil
      (cons nil (filter_sample (+ 1 count) nodes (cdr sample))))))

;Returns the Cooper and Herskovits value
;This is the function that should be called from code that uses this code
;    model should be described above
;    data should be in the format as described above
;    prior is the prior for the specified model. The equation specifies that
;           that the priors should be uniform.
(defun cooper (model data prior)
  (* (first_product (- (nodes data) 1) model data)
     prior))

;Returns the result of the left-most product in the formula
(defun first_product (nodes model data)
  (if (< nodes 0)
      1
    (* (second_product nodes (parent_instantiations nodes model data) data)
       (first_product (- nodes 1) model data))))

;Returns the result of the middle product in the formula
(defun second_product (node parent_instantiations data)
  (if (null parent_instantiations)
      1
    (/ (* (factorial (- (length (unique (states node data))) 1)) 
	  (third_product node (car parent_instantiations) (unique (states node data)) data)
	  (second_product node (cdr parent_instantiations) data))
       (factorial (- (+ (N (car parent_instantiations) data) (length (unique (states node data)))) 1)))))

;Returns the result of the right-most product in the formula 
(defun third_product (node parent_instantiation node_states data)
  (if (null node_states)
      1
    (* (factorial (alpha node (car node_states) parent_instantiation data))
       (third_product node parent_instantiation (cdr node_states) data))))