;;
;; FILE: learning.lsp
;; AUTH: Michael John Radwin <mjr@acm.org>
;;
;; DATE: Sat Nov  9 15:15:39 EST 1996
;;  $Id: learning.lsp,v 1.12 1997/03/27 20:49:18 mjr Exp mjr $
;;

(load (merge-pathnames (pathname "datapoint.lsp") *load-pathname*))
(load (merge-pathnames (pathname "fileio.lsp") *load-pathname*))

(defparameter *activation-epsillon*  0.01)
;(defparameter *activation-epsillon*  0.00001)
(defparameter *cluster-radius*       1)
;(defparameter *cluster-radius*       0.1)
(defparameter *activation-minimum*   0.1)
(defparameter *activation-prune*     0.9)
(defparameter *clusters*             nil)


(defun train-agglomerative (train-data stream)
  (format stream "~&;; train-agglomerative started at ~a (~a total cases)~%"
	  (get-universal-time) (length train-data))

  (setf *clusters* nil)
  (let ((i 0))
    (dolist (td train-data)
      (format t "~&;; introducing point ~a" i)
      (introduce-point td)
      (incf i)
      ))

  (format t "~&;; train-agglomerative finished at ~a~%"  (get-universal-time))
  )
  

(defmethod nn-point->cluster ((pt nn-point))
  (make-cluster :center (nn-point-point pt)
		:radius *cluster-radius*
		:class  (nn-point-class pt)))

;; we're assuming that this is a classification problem, so we'll
;; compare classes with eq instead of doing some sort of fancy
;; cross-validation.
(defmethod introduce-point ((pt nn-point))
  (let* ((class     (nn-point-class pt))
	 (activated 
	  (remove-if-not #'(lambda (e) (eq (cluster-class e) class))
			 (clusters-activated pt *activation-epsillon*))))
    (cond 
     ((null activated) 
;      (format t "~&;; none activated: new: ~a" (nn-point->cluster pt))
      (push (nn-point->cluster pt) *clusters*))
     ((= 1 (length activated))
;      (format t "~&;; one activated: ~a new: ~a" 
;	      (car activated) (nn-point->cluster pt))
      (setf *clusters* (remove (car activated) *clusters*))
      (push (merge-clusters (car activated) (nn-point->cluster pt)) *clusters*))
     (t
;      (format t "~&;; SOME ACTIVATED: ~a new: ~a" 
;	      activated (nn-point->cluster pt))
      (setf *clusters* (remove (car activated) *clusters*))
      (push (merge-clusters (car activated) (nn-point->cluster pt)) *clusters*)
      ))))


(defmethod agglomerative-predict ((unknown nn-point))
  ;; compute the activation
  (dolist (e *clusters*)
    (setf (cluster-activation e)
	  (exp (* -1 (distance (nn-point-point unknown)
			       (cluster-center e))))))

  ;; let's hope it's right!
  (let ((prediction (car (sort *clusters* #'> :key #'cluster-activation))))
    (values (cluster-class prediction) prediction)))


(defmethod clusters-activated ((unknown nn-point) epsillon)
  ;; we use the guassian for activation == exp(-d^2)
  (dolist (e *clusters*)
    (setf (cluster-activation e)
	  (exp (* -1 (distance (nn-point-point unknown)
			       (cluster-center e)))))
;    (format t "~&    CENTER: ~a ACTIVATION: ~a"
;	    (cluster-center e) (cluster-activation e))
    )
  (sort 
   (remove-if-not #'(lambda (e)
		      (>= (cluster-activation e)
			  epsillon))
		  *clusters*)
   #'> :key #'cluster-activation))

(defun inside-p (center radius query-point)
  (<= (distance center query-point) radius))

(defmethod merge-clusters ((e1 cluster) (e2 cluster))
  (format t "~&;; MERGING CLUSTERS: ~a~&;;                 : ~a" 
	  (list 'r (cluster-radius e1) 'ctr (cluster-center e1))
	  (list 'r (cluster-radius e2) 'ctr (cluster-center e2)))
;  (defun scale-point (coords scale)
;    (mapcar #'(lambda (i) (* i scale)) coords))

  (assert (eq (cluster-class e1) (cluster-class e2)))
  (if (inside-p (cluster-center e1) (cluster-radius e1) (cluster-center e2))
      (progn
	(format t "~&;;     INSIDE-P    : ~a"	
		(list 'r (cluster-radius e1) 'ctr (cluster-center e1)))
	e1)
    (let* ((vector (mapcar #'- (cluster-center e2) (cluster-center e1)))
	   (scale (/ (cluster-radius e2) (+ (cluster-radius e1) (cluster-radius e2))))
	   (scaled-v (mapcar #'(lambda (i) (* i scale)) vector))
	   (new-center (mapcar #'+ (cluster-center e1) scaled-v))
	   (new-cluster
	    (make-cluster :center new-center
			  :class (cluster-class e1)
			  :radius (max (cluster-radius e1)
				       (+ (distance new-center (cluster-center e2))
					  (cluster-radius e2))))))
      
;      (format t "~&;; scale      ~a" scale)
;      (format t "~&;; vector     ~a" vector)
;      (format t "~&;; scaled-v   ~a" scaled-v)
;      (format t "~&;; new-center ~a" new-center)
      
      (format t "~&;;     NEW CLUSTER : ~a"
	      (list 'r (cluster-radius new-cluster) 'ctr (cluster-center new-cluster)))
      
    ;; now send it back!
    new-cluster)))


;; take a raw datalist and split it into training and test data.  does
;; not modify nn-point-list.
(defun nn-point-list->train-and-test (nn-point-list &optional (test-pct .1))
  ;; make a shallow copy of the whole dataset, then remove random
  ;; elements, inserting them into test-data
  (let* ((train-data (mapcar #'identity nn-point-list))
	 (test-data '()))

    (setf *random-state* (make-random-state t))
    (dotimes (i (truncate (* test-pct (length train-data))))
      (let* ((num (random (length train-data)))
	     (elt (nth num train-data)))
;	(format t "~&i: ~a, len: ~a, rand: ~a, elt: ~a" 
;		i (length train-data) num elt)
	(push elt test-data)
	(delete elt train-data)))

    ;; two disjoint lists
    (values train-data test-data)))


;;------ distance and average functions for NN and k-NN ----------------

;; the distance between two d-dimensional points -- we use the square of
;; the geometric distance.
(defun delta-squared (x1 x2)
  (let ((delta (- x1 x2)))
    (* delta delta)))


(defun distance (point1 point2)
;  (assert (= (length point1) (length point2)))
  (reduce #'+
    (mapcar #'delta-squared point1 point2)))


(defmethod nn-point-distance ((p1 nn-point) (p2 nn-point))
  (distance (nn-point-point p1)
	    (nn-point-point p2)))


;; arithmetic mean.  plain 'n' simple.
(defun mean (&rest numbers)
  (/ (reduce #'+ numbers) 
     (length numbers)))


;; counts up symbols and returns the 'average' if there is one, or
;; :inconclusive if an average cannot be found.
(defun symbolic-mean (&rest symbols)
  ;; build an associative list of symbols and the times they appeared
  (let ((a-list '())
	(length (length symbols)))
    (dolist (sym symbols)
      (if (assoc sym a-list)
	  (incf (cdr (assoc sym a-list)))
	(setf a-list (acons sym 1 a-list))))

    ;; the mean is a symbol that appeared more than half the time
    (dolist (sym symbols)
      (when (> (/ (cdr (assoc sym a-list))
		  length)
	       0.5)
	(return-from symbolic-mean (car (assoc sym a-list)))))
    :inconclusive))



;;------ NN and k-NN API ----------------------------------------------

(defun make-unknown-nn-point (&rest coords)
  (make-nn-point
   :class nil
   :point coords))


;; vanilla nearest-neighbor.  returns 3 values:
;; 1. the class (or value) of the best point
;; 2. the best point itself
;; 3. the distance the best point was from the unknown point
(defmethod nearest-neighbor ((unknown nn-point) nn-point-list)
;  (format t "finding NN of sample: <~a, ~a> -- "
;	  (nn-point-point unknown) (nn-point-class unknown))
  (multiple-value-bind (neighbor point distance)
      (naive-nearest-neighbor unknown nn-point-list)
;    (format t "NN: <~a, ~a>" (nn-point-point point) neighbor)
    (values neighbor point distance)))
    
    


;; slightly better: k-nearest-neighbors.  like nearest-neighbor, it
;; returns 3 values, although the second and third are now lists.  The
;; first is the mean (or symbolic-mean) of the class/value of the best k
;; neighbors, which might possibly be :inconclusive
(defmethod k-nearest-neighbors ((unknown nn-point) nn-point-list k numeric?)
  (assert (< k (length nn-point-list)))
  (let ((neighbors 
	 (naive-k-nearest-neighbors unknown nn-point-list k)))
    (if numeric?
	(values (apply #'mean (mapcar #'nn-point-class neighbors))
		neighbors
		(mapcar #'nn-point-dist neighbors))
      (values (apply #'symbolic-mean (mapcar #'nn-point-class neighbors))
	      neighbors
	      (mapcar #'nn-point-dist neighbors)))))

(defmethod k-1357-nearest-neighbors ((unknown nn-point) nn-point-list numeric?)
  (let ((neighbors 
	 (naive-k-nearest-neighbors unknown nn-point-list 7)))
;    (format t "~&~a~%~a~%" unknown neighbors)
    (if numeric?
	(values 
	 (apply #'mean (mapcar #'nn-point-class         neighbors))
	 (apply #'mean (mapcar #'nn-point-class (subseq neighbors 0 5)))
	 (apply #'mean (mapcar #'nn-point-class (subseq neighbors 0 3)))
	 (nn-point-class                        (car    neighbors)))
      (values 
       (apply #'symbolic-mean (mapcar #'nn-point-class         neighbors))
       (apply #'symbolic-mean (mapcar #'nn-point-class (subseq neighbors 0 5)))
       (apply #'symbolic-mean (mapcar #'nn-point-class (subseq neighbors 0 3)))
       (nn-point-class                                 (car    neighbors))))))
      
  

;;------ naive nearest-neighbor algorithm -----------------------------

;; An O(dN) algorithm to find the nearest neighbor.  Returns 3 values:
;; 1. the class of the best point
;; 2. the best point itself
;; 3. the distance the best point was from the unknown point
(defmethod naive-nearest-neighbor ((unknown nn-point) nn-point-list)
  (let ((best-distance 9999999999999999)
	(best-point nil))
    (dolist (pt nn-point-list)
      (let ((dist (nn-point-distance unknown pt)))
;	(format t "~&dist: ~a, pt: ~a" dist (nn-point-point pt))
	(when (< dist best-distance)
	  (setf best-distance dist)
	  (setf best-point pt))))
    (values (nn-point-class best-point)
	    best-point
	    best-distance)))



;;------ naive k-nearest-neighbor algorithm ---------------------------

;; An O(N^2) algorithm to find the k nearest neighbors.  Returns a list
;; of nn-points of length k.
(defmethod naive-k-nearest-neighbors ((unknown nn-point) nn-point-list k)
  ;; side-effects nn-point-list
  (dolist (pt nn-point-list)
    (setf (nn-point-dist pt)
	  (nn-point-distance unknown pt)))

  ;; retrun the best k elements
  (subseq (sort nn-point-list #'< :key #'nn-point-dist) 0 k))