Matrix Example

extensions [ matrix ]

to setup
  clear-all

  let m1 matrix:from-column-list [[1 4 7][2 5 8][3 6 9]] ; 3x3 matrix
  let m2 matrix:from-row-list [[0 0 0][1 1 1][2 2 2][3 3 3]] ; 4x3 matrix
  let m3 matrix:make-identity 2 ; 2x2 identity matrix
  let m4 matrix:make-constant 5 1 100 ; 5x1 matrix, containing 100 in each entry

  output-print "m1"
  output-print m1
  output-print "m2"
  output-print m2
  output-print "m3"
  output-print m3
  output-print "m4"
  output-print m4

  output-print matrix:get m1 1 2   ; =>  (row 1, column 2), result is 6

  let mRef m2              ; mRef is another reference to the same matrix as m2
  let mCopy matrix:copy m2 ; => mCopy = new copy, not a reference to m2!

  matrix:set m2 0 1 9      ; Note mCopy remains unchanged, but m2 and mRef change

  output-print "m2 after matrix:set m2 0 1 9"
  output-print m2 ; => matrix: [[0 9 0][1 1 1][2 2 2][3 3 3]]
  output-print "mRef after matrix:set m2 0 1 9"
  output-print mRef ; same as above
  output-print "mCopy after matrix:set m2 0 1 9"
  output-print mCopy ; different -- no 9 in the first row.

  output-print "the dimensions of m4"
  output-print matrix:dimensions m4 ; => [5,1]
  output-print "m1 to a row list"
  output-print matrix:to-row-list m1 ; => [[1 2 3] [4 5 6] [7 8 9]]
  output-print "m4 to a column list"
  output-print matrix:to-column-list m4 ; => [[100][100][100][100][100]]

  let m5 matrix:plus-scalar m3 1
  output-print "m5 = m3 + scalar 1"
  output-print m5 ; => matrix:  [[2 1][1 2]]

  output-print "m3 + m5"
  output-print matrix:plus m3 m5 ; => [[3 1][1 3]]

  output-print "m5 x scalar 10"
  output-print matrix:times-scalar m5 10 ; => [[20 10][10 20]]
  output-print "matrix:times [[1 2][3 4]] [[0 1][-1 0]]"
  output-print matrix:times (matrix:from-row-list [[1 2][3 4]])
    (matrix:from-row-list [[0 1][-1 0]])  ; => {{matrix:  [ [ -2 1 ][ -4 3 ] ]}}

  output-print "matrix:times-element-wise [[1 2][3 4]] [[0 1][-1 0]]"
  output-print matrix:times-element-wise (matrix:from-row-list [[1 2][3 4]])
    (matrix:from-row-list [[0 1][-1 0]])  ; => {{matrix:    [[0 2][-3 0] ]}}

  output-print "matrix:inverse [[2 2][2 0]]"
  output-print matrix:inverse (matrix:from-row-list [[2 2][2 0]])  ; => {{matrix:  [ [ 0 0.5 ][ 0.5 -0.5 ] ]}}
  carefully [
    output-print matrix:inverse matrix:from-row-list [[0 0] [0 0]]
  ][
    output-print "Can't invert matrix: [[0 0] [0 0]]!"
  ]

  output-print "matrix:transpose m2"
  output-print matrix:transpose m2 ; => matrix: [[0 1 2 3][9 1 2 3][0 1 2 3]]
  output-print "matrix:submatrix m1 0 1 2 3"
  output-print matrix:submatrix m1 0 1 2 3 ; matrix, rowStart, colStart, rowEnd, colEnd
                                           ; rows from 0 (inclusive) to 2 (exclusive),
                                           ; columns from 1 (inclusive) to 3 (exclusive)
                                           ; => matrix: [[2 3][5 6]]
  output-print "matrix:get-row m1 2"
  output-print matrix:get-row m1 1         ; puts the second row of m1 (starting from 0) in a simple
                                           ; NetLogo list.  => [4 5 6]
  output-print "matrix:get-column m1 2"
  output-print matrix:get-column m1 1      ; puts the second column of m1 (starting from 0) in a simple
                                           ; NetLogo list.  => [2 5 8]

  output-print "create a new matrix m6 from row list [[0 0 0][1 1 1][2 2 2][3 3 3]]"
  let m6 matrix:from-row-list [[0 0 0][1 1 1][2 2 2][3 3 3]] ; 4x3 matrix
  output-print m6
  output-print "matrix:set-row m6 1 [20 21 22]"
  matrix:set-row m6 1 [20 21 22]  ; => [[0 0 0][20 21 22][2 2 2][3 3 3]]
  output-print m6
  output-print "matrix:set-column m6 2 [10 11 12 13]"
  matrix:set-column m6 2 [10 11 12 13]  ; => [[0 0 10][20 21 11][2 2 12][3 3 13]]
  output-print m6
  output-print "matrix:swap-rows m6 1 2"
  matrix:swap-rows m6 1 2  ; => [[0 0 10][2 2 12][20 21 11][3 3 13]]
  output-print m6
  output-print "matrix:swap-columns m6 2 0"
  matrix:swap-columns m6 2 0  ; => [[10 0 0][12 2 2][11 21 20][13 3 3]]
  output-print m6

  output-print "matrix:real-eigenvalues m5"
  ; recall m5 = [[2 1][1 2]]
  output-print matrix:real-eigenvalues m5 ; => [1.0000000000000002 3] NOTE: accuracy not perfect.
  output-print "matrix:imaginary-eigenvalues m5"
  output-print matrix:imaginary-eigenvalues m5 ; => [ 0 0 ]
  output-print "matrix:eigenvectors m5"
  output-print matrix:eigenvectors m5 ; => matrix:  [[0.707... 0.707...][-0.707... 0.707...]] (.707... = sqrt(2)/2 )

  let A matrix:from-row-list [[1 2][3 4]]
  let B matrix:from-row-list [[-2 1][-4 3]]
  output-print "matrix:solve [[1 2][3 4]] [[-2 1][-4 3]]"
  output-print matrix:solve A B ; => [[0 1][-1 0]]
                                ;  Solve is sort of like matrix division.
                                ;  That is, it tries to find matrix X such that  A * X = B.
                                ;  It gives a "least-squares" solution, if no perfect solution exists.

  output-print "matrix:set-and-report m1 2 2 0"
  output-print matrix:set-and-report m1 2 2 0  ;; m1 is unchanged, a new matrix is reported.
  output-print "m1 is unchanged"
  output-print m1
  set m1 matrix:set-and-report m1 2 2 0        ;; can also change m1 itself, using the reporter.

  output-print "matrix:forecast-linear-growth [20 25 28 32 35 39]  => ?? linear growth rate"
  output-print matrix:forecast-linear-growth [20 25 28 32 35 39]  ; a linear extrapolation of the next item in the list.

  output-print "matrix:forecast-compound-growth [20 25 28 32 35 39] => ~13.7% compound growth rate"
  output-print matrix:forecast-compound-growth [20 25 28 32 35 39]
  output-print "matrix:forecast-compound-growth [2 2.2 2.42 2.662 2.9282] => 10.0% compound growth rate"
  output-print matrix:forecast-compound-growth [2 2.2 2.42 2.662 2.9282]
  output-print "matrix:forecast-compound-growth [4 3.6 3.24 2.916 2.6244] => minus 10.0% compound growth rate"
  output-print matrix:forecast-compound-growth [4 3.6 3.24 2.916 2.6244]

  output-print "matrix:forecast-continuous-growth [20 25 28 32 35 39] => ~12.8% continuous growth rate"
  output-print matrix:forecast-continuous-growth [20 25 28 32 35 39]
  output-print "matrix:forecast-continuous-growth [2 2.2 2.42 2.662 2.9282] => ~9.5% continuous growth rate"
  output-print matrix:forecast-continuous-growth [2 2.2 2.42 2.662 2.9282] ; a compound growth extrapolation of the next item in the list.
  output-print "matrix:forecast-continuous-growth [4 3.6 3.24 2.916 2.6244] => ~ minus 10.5% continuous growth rate"
  output-print matrix:forecast-continuous-growth [4 3.6 3.24 2.916 2.6244]

  output-print "matrix:regress matrix:from-column-list [[2 4 5 8 10] [3 4 3 7 8] [2 3 5 8 9]]"
  output-print matrix:regress matrix:from-column-list [[2 4 5 8 10] [3 4 3 7 8] [2 3 5 8 9]]
  output-print "matrix:regress matrix:from-row-list [[2 3 2] [4 4 3] [5 3 5] [8 7 8] [10 8 9]]"
  output-print matrix:regress matrix:from-row-list [[2 3 2] [4 4 3] [5 3 5] [8 7 8] [10 8 9]]

  ; The forecasts may also be done using regress.  Here they are. Note that the next value
  ; to be forecast is for the value after the 6 observations supplied. But since lists begin with position
  ; zero, that will be for position 6.
  ; this is what happens inside matrix:forecast-linear-growth.
  output-print "long way to do the linear forecast using regress"
  let data-list [20 25 28 32 35 39]
  let indep-var (n-values length data-list [?]) ; 0,1,2...,5
  let lin-output matrix:regress matrix:from-column-list (list data-list indep-var)
  let lincnst item 0 (item 0 lin-output)
  let linslpe item 1 (item 0 lin-output)
  let linR2   item 0 (item 1 lin-output)
  output-print (list (lincnst + linslpe * 6) (lincnst) (linslpe) (linR2))

  ; this is what happens inside matrix:forecast-compound-growth.
  output-print "long way to do the compound-growth forecast using regress"
  let com-log-data-list  (map [ln ?] [20 25 28 32 35 39])
  let com-indep-var2 (n-values length com-log-data-list [?]) ; 0,1,2...,5
  let com-output matrix:regress matrix:from-column-list (list com-log-data-list com-indep-var2)
  let comcnst exp item 0 (item 0 com-output)
  let comrate exp item 1 (item 0 com-output)
  let comR2       item 0 (item 1 com-output)
  output-print (list (comcnst * comrate ^ 6) (comcnst) (comrate) (comR2))

  ; this is what happens inside matrix:forecast-continuous-growth.
  output-print "long way to do the continuous-growth forecast using regress"
  let con-log-data-list  (map [ln ?] [20 25 28 32 35 39])
  let con-indep-var2 (n-values length con-log-data-list [?]) ; 0,1,2...,5
  let con-output matrix:regress matrix:from-column-list (list con-log-data-list con-indep-var2)
  let concnst exp item 0 (item 0 con-output)
  let conrate     item 1 (item 0 con-output)
  let conR2       item 0 (item 1 con-output)
  output-print (list (concnst * exp (conrate * 6)) (concnst) (conrate) (conR2))


  ;;NOTE: The matrix extension could be extended with more functionality,
  ;;      since the Jama library it is based on can do much more
  ;;        (e.g. LU, Cholesky, SV decompositions)
end 


; Public Domain:
; To the extent possible under law, Uri Wilensky has waived all
; copyright and related or neighboring rights to this model.