source: project/release/4/signal-diagram/trunk/examples/Izhikevich03.scm @ 30942

;;
;; Simple Model of Spiking Neurons, Izhikevich  E
;;
;; IEEE Transactions on Neural Networks (2003) 14:1569- 1572
;; 
;; 

(use signal-diagram signal-diagram-dynamics mathh srfi-1 )

(define (Izhikevich03:construct t V U k1 k2 k3 theta a b c d spike tspike Isyn )

  `((f . ,(let* ((subthreshold-eq         (ODE (variable h) t 
                                               (V  (+ (* k1 V V)  (* k2 V) k3 (neg U) Isyn))
                                               (U  (* a (- (* b V) U)))))

                 (threshold-detect        (ASSIGN (spike (> V theta))))

                 (tspike-set              (ASSIGN (tspike (if spike t tspike))))

                 (subthreshold-regime     (SEQUENCE subthreshold-eq (SEQUENCE threshold-detect tspike-set)))
                 
                 (refractory-eq           (ASSIGN (t (+ t h) )
                                                  (spike false)
                                                  (tspike tspike)
                                                  (V c)
                                                  (U (+ U d))))
                 )

            (TRANSIENT subthreshold-regime refractory-eq 'spike)

            ))

    (initial .
             ((h              1e-2)
             (t               ,t)
             (V               ,V)
             (U               ,U)
             (k1              ,k1)
             (k2              ,k2)
             (k3              ,k3)
             (theta           ,theta)
             (a               ,a)
             (b               ,b)
             (c               ,c)
             (d               ,d)
             (spike           ,spike)
             (tspike          ,tspike)
             (Isyn            ,Isyn)
             ))
    ))


(define (dataflow model)
  (let* ((f         (alist-ref 'f model)) 
         (initial   (alist-ref 'initial model))
         (input     (map car initial)))
    (print "input = " input)
    (print "dataflow f: "  (dataflow (construct f) input))
    ))


(define (codegen name model #!key (language 'scheme) (random #f) (pre #t) (post #t) (solver 'rkdp))
  (let* ((f         (alist-ref 'f model)) 
         (initial   (alist-ref 'initial model)))
    (case language
      ((octave Octave) (codegen/Octave name (construct f) initial: initial pre: pre solver: solver))
      ((scheme Scheme) (codegen/scheme name (construct f) initial: initial pre: pre solver: solver))
      ((ML)     (codegen/ML name (construct f) initial: initial random: random pre: pre post: post solver: solver)))
    ))


;; Parameters common to all spiking regimes
(define k1    0.04)
(define k2    5.0)
(define k3    140.0)
(define theta 30.0)
(define Isyn  10.0)

;; State initial values
(define t 0.0)
(define V -65.0)


;; Regular spiking (RS) parameters
(define RS 
  (let ((a 0.02)
        (b 0.2)
        (c -65.0)
        (d 8.0))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))


;; Intrinsically bursting (IB) parameters
(define IB 
  (let ((a 0.02)
        (b 0.2)
        (c -55.0)
        (d 4.0))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))


;; Chattering (CH) parameters
(define CH
  (let ((a 0.02)
        (b 0.2)
        (c -50.0)
        (d 2.0))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))


;; Fast-spiking (FS) parameters
(define FS
  (let ((a 0.1)
        (b 0.2)
        (c -65.0)
        (d 2.0))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))



;; Thalamo-cortical (TC) parameters
(define TC
  (let ((a 0.02)
        (b 0.25)
        (c -65.0)
        (d 0.05))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))


;; Resonator (RZ) parameters
(define RZ
  (let ((a 0.1)
        (b 0.26)
        (c -65.0)
        (d 2.0))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))



;; Low-threshold spiking (LTS) parameters
(define LTS
  (let ((a 0.02)
        (b 0.25)
        (c -65.0)
        (d 2.0))
    (Izhikevich03:construct
     ;; t V U k1 k2 k3 theta a b c d spike tspike Isyn 
     t V (* b V) k1 k2 k3 theta a b c d #f t Isyn)))


(define models `((RS  . ,RS) (IB . ,IB) (CH . ,CH) (FS . ,FS)  (RZ . ,RZ) (LTS . ,LTS) ))

  
(with-output-to-file "Izhikevich03_solver.sml"
  (lambda () 
    (let recur ((models models) (pre #t))
      (if (pair? models)
          (let* ((name.model (car models))
                 (name  (car name.model))
                 (model (cdr name.model)))
            (codegen name model language: 'ML solver: 'rkdp random: #f pre: pre post: (null? (cdr models)))
            (recur (cdr models) #f)
            )))
    ))

(with-output-to-file "Izhikevich03_solver.m"
  (lambda () 
    (let recur ((models models) (pre #t))
      (if (pair? models)
          (let* ((name.model (car models))
                 (name  (car name.model))
                 (model (cdr name.model)))
            (codegen name model language: 'octave pre: pre solver: 'lsode)
            (recur (cdr models) #f)
            )))
    ))

(with-output-to-file "Izhikevich03_solver.scm"
  (lambda () 
    (let recur ((models models) (pre #t))
      (if (pair? models)
          (let* ((name.model (car models))
                 (name  (car name.model))
                 (model (cdr name.model)))
            (codegen name model language: 'scheme pre: pre)
            (recur (cdr models) #f)
            )))
    ))