logistic, ridge classifier

logistic.clj

(ns forma.source.logistic
  (:use [clojure.math.numeric-tower :only (sqrt floor abs expt)]
        [clojure-csv.core])
  (:require [incanter.core :as i]))

(defn feature-vec
  [n]
  (map (partial cons 1)
       (for [x (range n)]
         (take 22 (repeatedly rand)))))

(defn label-vec
  [n]
  (for [x (range n)] (if (> (rand) 0.5) 1 0)))

(defn scaled-vector
  [scalar coll]
  (map #(* scalar %) coll))

(defn logistic-fn
  [x]
  (let [exp-x (Math/exp x)]
    (/ exp-x (inc exp-x))))

(defn dot-product
  [x y]
  (reduce + (map * x y)))

(defn logistic-prob
  [beta-seq x]
  (logistic-fn (dot-product beta-seq x)))

(defn log-likelihood
  [beta-seq label x]
  (let [prob (logistic-prob beta-seq x)]
    (+ (* label (Math/log prob))
       (* (- 1 label) (Math/log (- 1 prob))))))

(defn total-log-likelihood
  "returns the total log likelihood for a group of pixels; input
  labels and features for the group of pixels, aligned correctly so
  that the first label and feature correspond to the first pixel."
  [beta-seq labels features]
  (reduce + (map (partial log-likelihood beta-seq) labels features)))

(defn probability-calc
  "returns a vector of probabilities for each observation"
  [beta-seq feature-seqs]
  (map (partial logistic-prob beta-seq)
       feature-seqs))

(defn score-seq
  "returns the scores for each parameter"
  [beta-seq labels features]
  (let [prob-seq (probability-calc beta-seq features)]
    (i/mmult (i/trans features) (map - labels prob-seq))))

(defn info-matrix
  "returns the information matrix for the logistic probability
  function"
  [beta-seq labels features]
  (let [mult-func (fn [x] (* x (- 1 x)))
        prob-seq  (->> (probability-calc beta-seq features)
                       (map mult-func))]
    (i/mmult (map scaled-vector prob-seq (i/trans features))
             features)))

(defn beta-update
  [beta-seq labels features rdg-cons]
  (let [num-features (count beta-seq)
        info-adj (i/plus
                  (info-matrix beta-seq labels features)
                  (i/diag (repeat num-features rdg-cons)))]
    (i/mmult
     (i/solve info-adj)
     (score-seq beta-seq labels features))))

(defn logistic-beta-vector 
  [labels features rdg-cons]
  (let [b (repeat 23 0)]
    (loop [beta b iter 9]
      (if (zero? iter)
        beta
        (recur (let [inc-beta (beta-update beta labels features rdg-cons)]
                 (println "")
                 (println (flatten beta))
                 map + beta inc-beta)
               (dec iter))))))


;; (time (dotimes [x 1]
;;         (doall
;;          (logistic-beta-vector
;;           (label-vec 100000)
;;           (feature-vec 100000)
;;           0.00000001))))