annotate yeti/silvet.yeti @ 20:982aa1197a7e

Getting there, slowly, sort of, with EM
author Chris Cannam
date Thu, 27 Mar 2014 11:52:07 +0000
parents f1f8c84339d0
children 8e61ec97b34e
rev   line source
Chris@12 1
Chris@12 2 program silvet;
Chris@12 3
Chris@12 4 { prepareTimeFrequency } = load timefreq;
Chris@13 5 { loadTemplates, extractRanges } = load templates;
Chris@14 6
Chris@14 7 em = load em;
Chris@19 8 em1 = load em_onecolumn;
Chris@12 9
Chris@12 10 mat = load may.matrix;
Chris@12 11 vec = load may.vector;
Chris@14 12 plot = load may.plot;
Chris@12 13
Chris@12 14 templates = loadTemplates ();
Chris@12 15
Chris@13 16 ranges = extractRanges templates;
Chris@13 17
Chris@13 18 eprintln "\nWe have \(length (keys templates)) instruments:";
Chris@13 19 for (sort (keys templates)) do k:
Chris@14 20 eprintln " * \(k) \(mat.size templates[k]) range \(ranges[k].lowest) -> \(ranges[k].highest)";
Chris@13 21 done;
Chris@12 22 eprintln "";
Chris@12 23
Chris@12 24 columns = prepareTimeFrequency "test.wav";
Chris@12 25
Chris@19 26 /*
Chris@12 27 height = if empty? columns then 0 else vec.length (head columns) fi;
Chris@12 28
Chris@13 29 chunkSize = { rows = height, columns = 100 };
Chris@13 30
Chris@14 31 emdata = em.initialise ranges templates 88 chunkSize;
Chris@13 32
Chris@14 33 eprintln "initialised EM data: overall pitch range \(emdata.lowest) -> \(emdata.highest)";
Chris@13 34
Chris@12 35 chunkify cols =
Chris@12 36 if empty? cols then []
Chris@12 37 else
Chris@13 38 (mat.resizedTo chunkSize
Chris@13 39 (mat.fromColumns (take chunkSize.columns cols)))
Chris@13 40 :. \(chunkify (drop chunkSize.columns cols));
Chris@12 41 fi;
Chris@12 42
Chris@12 43 chunks = chunkify columns;
Chris@12 44
Chris@12 45 eprintln "we have \(length chunks) chunks of size \(mat.size (head chunks))";
Chris@12 46
Chris@14 47 eprintln "attempting one expectation phase...";
Chris@13 48
Chris@14 49 error = em.performExpectation emdata (head chunks);
Chris@14 50
Chris@14 51 eprintln "done, result has dimension \(mat.size error)";
Chris@14 52
Chris@14 53 eprintln "attempting one maximisation phase...";
Chris@14 54
Chris@14 55 newP = em.performMaximisation emdata (head chunks) error;
Chris@14 56
Chris@14 57 eprintln "done";
Chris@14 58
Chris@14 59 \() (plot.plot [ Grid (head chunks) ]);
Chris@14 60 \() (plot.plot [ Grid error ]);
Chris@14 61
Chris@14 62 \() (plot.plot [ Grid newP ]);
Chris@19 63 */
Chris@19 64
Chris@19 65 em1data = em1.initialise ranges templates 88;
Chris@19 66
Chris@19 67 col = head (drop 35 columns);
Chris@19 68
Chris@20 69 \() (plot.plot [ Caption "Source frequency distribution", Vector col ]);
Chris@19 70
Chris@20 71 sourceGrid d =
Chris@20 72 mat.fromColumns (map do k: d.sources[k] done (sort (keys d.sources)));
Chris@20 73
Chris@20 74 \() (plot.plot [ Caption "Source distribution beforehand", Grid (sourceGrid em1data)]);
Chris@20 75
Chris@20 76 oneIteration em1data col n =
Chris@19 77 (q = em1.performExpectation em1data col;
Chris@20 78 // \() (plot.plot [ Caption "Frequency distribution and output of E-step for iteration \(n)", Vector col, Vector q ]);
Chris@19 79 newdata = em1.performMaximisation em1data col q;
Chris@20 80 if (n % 2 == 0) then
Chris@20 81 \() (plot.plot [ Caption "Pitch distribution before and after M-step update for iteration \(n)", Vector (em1data.pitches), Vector (newdata.pitches) ]);
Chris@20 82 \() (plot.plot [ Caption "Source distribution after M-step update for iteration \(n)", Grid (sourceGrid newdata) ]);
Chris@20 83 fi;
Chris@19 84 newdata);
Chris@19 85
Chris@20 86 iterations = 6;
Chris@19 87
Chris@20 88 var d = em1data;
Chris@20 89
Chris@20 90 for [1..iterations] do i:
Chris@20 91 d := oneIteration d col i;
Chris@20 92 done;
Chris@14 93
Chris@14 94 ();
Chris@14 95
Chris@14 96