GRN (gene regulatory network) experimental discrete formulation

core.clj

(ns org.nfrac.gron.core
  (:require [clojure.spec :as s]
            [clojure.spec.gen :as gen]))

(s/def ::sep #{-1})

(s/def ::base (-> nat-int?
                  (s/with-gen #(s/gen (s/int-in 0 10)))))

(s/def ::seq-gene
  (s/cat :sep1 ::sep
         ::promoter (s/* ::base)
         :sep2 ::sep
         ::repressor (s/* ::base)
         :sep3 ::sep
         ::produce (s/* ::base)
         :sep4 ::sep
         ::consume (s/* ::base)))

(s/def ::genome
  (s/+ ::seq-gene))

(s/def ::grn
  (s/keys :req [::genome]))

(comment
 (s/exercise ::seq-gene))

(s/def ::base-counts
  (s/map-of ::base (-> nat-int? (s/with-gen #(s/gen (s/int-in 1 4))))))

(s/def ::cell ::base-counts)

(s/def ::promoter ::base-counts)
(s/def ::repressor ::base-counts)
(s/def ::produce ::base-counts)
(s/def ::consume ::base-counts)
(s/def ::gene (s/keys :req [::promoter
                            ::repressor
                            ::produce
                            ::consume]))

(comment
 (s/exercise ::gene))

(defn regulator-matched?
  [regulator cell]
  (and (seq regulator)
       (every? (fn [[base n]] (>= (cell base 0) n)) regulator)))

(s/fdef regulator-matched?
        :args (s/cat :regulator ::base-counts
                     :cell ::cell))

(defn gene-expressed?
  [gene cell]
  (and (regulator-matched? (::promoter gene) cell)
       (not (regulator-matched? (::repressor gene) cell))))

(s/fdef gene-expressed?
        :args (s/cat :gene ::gene
                     :cell ::cell))

(defn expressed-genes
  [genes cell]
  (loop [genes genes
         bound-cell cell
         expressed ()]
    (if-let [gene (first genes)]
      (if (gene-expressed? gene bound-cell)
        (recur (rest genes)
          (merge-with - bound-cell (::promoter gene) (::repressor gene))
          (conj expressed gene))
        ;; gene not expressed
        (recur (rest genes) bound-cell expressed))
      ;; checked all genes
      expressed)))

(s/fdef expressed-genes
        :args (s/cat :genes (s/every ::gene)
                     :cell ::cell)
        :ret (s/every ::gene))

(defn gene-effects
  [genes]
  (let [pro (apply merge-with + (map ::produce genes))
        con (apply merge-with + (map ::consume genes))]
    (merge-with + pro (zipmap (keys con) (map - (vals con))))))

(s/fdef gene-effects
        :args (s/cat :genes (s/every ::gene))
        :ret (s/every-kv ::base int?))

(defn step
  [genes cell inputs]
  (let [icell (merge-with + cell inputs)
        xg (expressed-genes genes icell)
        ge (gene-effects xg)]
    (reduce (fn [cell [base delta]]
              (assoc cell base
                     (max 0 (+ (cell base 0) delta))))
            cell ;icell
            ge)))

(s/fdef step
        :args (s/cat :genes (s/every ::gene)
                     :cell ::cell
                     :inputs ::base-counts)
        :ret ::cell)

(defn genes-form
  [seq-genome]
  (->> seq-genome
       (s/conform ::genome)
       (map (fn [m]
              (reduce (fn [m k]
                        (assoc m k (frequencies (get m k))))
                      m
                      [::promoter ::repressor ::produce ::consume])))))

(s/fdef genes-form
        :args (s/cat :seq-genome ::genome)
        :ret (s/every ::gene))

(comment
 (require 'proto-repl-charts.charts)
 (let [genes (genes-form (gen/generate (s/gen ::genome)))
       cell-ts (reductions (fn [cell in]
                             (step genes cell in))
                           (gen/generate (s/gen ::base-counts))
                           (repeat 30 {}))
       bases (range 10)
       base-ts-map (reduce (fn [m b]
                             (assoc m b (map #(get % b 0) cell-ts)))
                           {}
                           bases)]
   (println (count genes) "genes")

   (proto-repl-charts.charts/line-chart
     "Random GRNs"
     base-ts-map
     {:labels (range (count cell-ts))})))