001 package com.croftsoft.apps.cyborg;
002
003 import java.io.*;
004 import java.util.*;
005
006 import org.apache.commons.math.MathException;
007 import org.apache.commons.math.distribution.DistributionFactory;
008 import org.apache.commons.math.distribution.NormalDistribution;
009 import org.apache.commons.math.stat.descriptive.DescriptiveStatistics;
010
011 /*********************************************************************
012 * Surrogate directional activity.
013 *
014 * Translated from SurrogateDirectionalActivityV2.xls (2004 October)
015 * by Lawrence J. Cauller, Ph.D., which was donated to the
016 * Public Domain by Dr. Cauller on 2004-12-01.
017 *
018 * @version
019 * 2004-12-23
020 * @since
021 * 2004-10-19
022 * @author
023 * <a href="https://www.croftsoft.com/">David Wallace Croft</a>
024 *********************************************************************/
025
026 public final class Activity
027 //////////////////////////////////////////////////////////////////////
028 //////////////////////////////////////////////////////////////////////
029 {
030
031 private static final int DIRECTIONS = 8;
032
033 private static final int NEURONS = 24;
034
035 private static final int TRIALS = 100;
036
037 //
038
039 private static final double RATE_SPONTANEOUS = 10.0;
040
041 private static final double RATE_MIN = 5.0;
042
043 private static final double RATE_MAX = 50.0;
044
045 private static final double RATE_DEVIATION = 10.0;
046
047 //
048
049 /** Random number generator. */
050 private static final Random random = new Random ( );
051
052 /** Used to create normal distributions for adding noise. */
053 private static final DistributionFactory distributionFactory
054 = DistributionFactory.newInstance ( );
055
056 //
057
058 private static final String OUTPUT_FILENAME = "out.txt";
059
060 //////////////////////////////////////////////////////////////////////
061 //////////////////////////////////////////////////////////////////////
062
063 public static void main ( String [ ] args )
064 throws Exception
065 //////////////////////////////////////////////////////////////////////
066 {
067 // Write the output to a file using PrintStream printf().
068
069 PrintStream printStream = new PrintStream (
070 new FileOutputStream ( OUTPUT_FILENAME ) );
071
072 // Assume we know the preferred stimulus direction for each of the
073 // neurons that we are recording from.
074
075 double [ ] preferredDirections = new double [ NEURONS ];
076
077 for ( int unit = 0; unit < NEURONS; unit++ )
078 {
079 // Assume that the neurons we chose have evenly distributed
080 // preferred stimulus directions
081
082 preferredDirections [ unit ]
083 = ( 2.0 * Math.PI * ( unit ) ) / NEURONS;
084 }
085
086 printTemplateActivations ( printStream, preferredDirections );
087
088 runTrials ( printStream, preferredDirections );
089 }
090
091 /*********************************************************************
092 * Prints the template activations with one standard dev of noise.
093 *********************************************************************/
094 public static void printTemplateActivations (
095 PrintStream printStream,
096 double [ ] preferredDirections )
097 throws Exception
098 //////////////////////////////////////////////////////////////////////
099 {
100 double [ ] templateActivations = new double [ NEURONS ];
101
102 double [ ] addDeviations = new double [ NEURONS ];
103
104 double [ ] subDeviations = new double [ NEURONS ];
105
106 printStream.println ( "unit preferred template addDev subDev" );
107
108 for ( int unit = 0; unit < NEURONS; unit++ )
109 {
110 templateActivations [ unit ] = computeNoiselessActivation (
111 0.0, preferredDirections [ unit ] );
112
113 addDeviations [ unit ]
114 = templateActivations [ unit ] + RATE_DEVIATION;
115
116 subDeviations [ unit ]
117 = templateActivations [ unit ] - RATE_DEVIATION;
118
119 subDeviations [ unit ] = Math.ceil ( subDeviations [ unit ] );
120
121 if ( subDeviations [ unit ] <= 0.0 )
122 {
123 subDeviations [ unit ] = 0.0;
124 }
125
126 printStream.printf (
127 "%1$02d %2$03.0f %3$06.3f %4$06.3f %5$06.3f%n",
128 new Object [ ] {
129 new Integer ( unit ),
130 new Double ( Math.toDegrees (
131 preferredDirections [ unit ] ) ),
132 new Double ( templateActivations [ unit ] ),
133 new Double ( addDeviations [ unit ] ),
134 new Double ( subDeviations [ unit ] ) } );
135 }
136 }
137
138 /*********************************************************************
139 * Runs the trials, gathers statistics, and prints the results.
140 *********************************************************************/
141 public static void runTrials (
142 PrintStream printStream,
143 double [ ] preferredDirections )
144 throws Exception
145 //////////////////////////////////////////////////////////////////////
146 {
147 // initialize statistics objects
148
149 DescriptiveStatistics [ ] descriptiveStatisticsArray
150 = new DescriptiveStatistics [ DIRECTIONS ];
151
152 DescriptiveStatistics [ ] [ ]
153 averagePopulationVectorsDescriptiveStatistics
154 = new DescriptiveStatistics [ DIRECTIONS ] [ NEURONS ];
155
156 for ( int i = 0; i < DIRECTIONS; i++ )
157 {
158 descriptiveStatisticsArray [ i ]
159 = DescriptiveStatistics.newInstance ( );
160
161 for ( int j = 0; j < NEURONS; j++ )
162 {
163 averagePopulationVectorsDescriptiveStatistics [ i ] [ j ]
164 = DescriptiveStatistics.newInstance ( );
165 }
166 }
167
168 DescriptiveStatistics grandDescriptiveStatistics
169 = DescriptiveStatistics.newInstance ( );
170
171 // run the trial data
172
173 printStream.println ( "\n" + "trial actual predicted error" );
174
175 for ( int trial = 0; trial < TRIALS; trial++ )
176 {
177 int directionIndex = ( trial * DIRECTIONS ) / TRIALS;
178
179 double actualDirection
180 = directionIndex * 2.0 * Math.PI / DIRECTIONS;
181
182 double [ ] noisyActivations = new double [ NEURONS ];
183
184 for ( int unit = 0; unit < NEURONS; unit++ )
185 {
186 // record the activity from a single neuron for the given trial
187
188 noisyActivations [ unit ] = generateNoisyActivation (
189 preferredDirections [ unit ],
190 actualDirection );
191
192 // record the statistical data for this direction and this unit
193
194 averagePopulationVectorsDescriptiveStatistics
195 [ directionIndex ] [ unit ].addValue (
196 noisyActivations [ unit ] );
197 }
198
199 // Knowing only the preferred directions for each neuron and their
200 // current activations, predict the direction for this trial.
201
202 double predictedDirection = predictDirection (
203 preferredDirections,
204 noisyActivations );
205
206 // Compute the error between predicted and actual direction
207
208 double error = Math.abs ( Math.min (
209 predictedDirection - actualDirection,
210 actualDirection + 2.0 * Math.PI - predictedDirection ) );
211
212 // How variable is the error for the trials in a given direction?
213
214 descriptiveStatisticsArray [ directionIndex ].addValue ( error );
215
216 // How variable is the error from trial to trial?
217
218 grandDescriptiveStatistics.addValue ( error );
219
220 // trial actual predicted error
221
222 printStream.printf (
223 "%1$02d %2$03.0f %3$07.3f %4$07.3f%n",
224 new Object [ ] {
225 new Integer ( trial ),
226 new Double ( Math.toDegrees ( actualDirection ) ),
227 new Double ( Math.toDegrees ( predictedDirection ) ),
228 new Double ( Math.toDegrees ( error ) ) } );
229 }
230
231 printStream.println ( "\n" +
232 "direction / mean error / error std / mean-error / mean+error" );
233
234 for ( int i = 0; i < DIRECTIONS; i++ )
235 {
236 double direction = i * 2.0 * Math.PI / DIRECTIONS;
237
238 double meanError = descriptiveStatisticsArray [ i ].getMean ( );
239
240 double errorStd
241 = descriptiveStatisticsArray [ i ].getStandardDeviation ( );
242
243 // the mean error and standard deviation in error for a direction
244
245 printStream.printf (
246 "%1$03.0f %2$06.3f %3$06.3f %4$06.3f %5$06.3f%n",
247 new Double ( Math.toDegrees ( direction ) ),
248 new Double ( Math.toDegrees ( meanError ) ),
249 new Double ( Math.toDegrees ( errorStd ) ),
250 new Double ( Math.toDegrees ( meanError - errorStd ) ),
251 new Double ( Math.toDegrees ( meanError + errorStd ) ) );
252 }
253
254 // the mean error and standard deviation in error for all trials
255
256 printStream.println ( "\n" + "grand mean error: "
257 + Math.toDegrees ( grandDescriptiveStatistics.getMean ( ) ) );
258
259 printStream.println ( "\n" + "grand deviation: "
260 + Math.toDegrees (
261 grandDescriptiveStatistics.getStandardDeviation ( ) ) );
262
263 // print the mean activity from each neuron for each direction
264
265 printStream.println ( "\n" + "direction preferred average" );
266
267 for ( int i = 0; i < DIRECTIONS; i++ )
268 {
269 for ( int j = 0; j < NEURONS; j++ )
270 {
271 printStream.printf (
272 "%1$03.0f %2$03.0f %3$03.0f%n",
273 new Double [ ] {
274 new Double (
275 Math.toDegrees ( i * 2.0 * Math.PI / DIRECTIONS ) ),
276 new Double (
277 Math.toDegrees ( j * 2.0 * Math.PI / NEURONS ) ),
278 new Double (
279 Math.ceil (
280 averagePopulationVectorsDescriptiveStatistics
281 [ i ] [ j ]
282 .getMean ( ) ) ) } );
283 }
284 }
285 }
286
287 /*********************************************************************
288 * Generates activation for a neuron without random noise.
289 *********************************************************************/
290 public static double computeNoiselessActivation (
291 double preferredDirection,
292 double actualDirection )
293 //////////////////////////////////////////////////////////////////////
294 {
295 double angleFromPreferred
296 = preferredDirection - actualDirection;
297
298 if ( angleFromPreferred < 0.0 )
299 {
300 angleFromPreferred += 2.0 * Math.PI;
301 }
302
303 double computedActivation;
304
305 if ( ( angleFromPreferred <= Math.PI / 2 )
306 || ( angleFromPreferred >= 3 * Math.PI / 2 ) )
307 {
308 computedActivation
309 = RATE_SPONTANEOUS
310 + Math.cos ( angleFromPreferred )
311 * ( RATE_MAX - RATE_SPONTANEOUS );
312 }
313 else
314 {
315 computedActivation
316 = RATE_SPONTANEOUS
317 + Math.cos ( angleFromPreferred )
318 * ( RATE_SPONTANEOUS - RATE_MIN );
319 }
320
321 return computedActivation;
322 }
323
324 /*********************************************************************
325 * Generates activation for a neuron with random noise.
326 *********************************************************************/
327 public static double generateNoisyActivation (
328 double preferredDirection,
329 double actualDirection )
330 throws Exception
331 //////////////////////////////////////////////////////////////////////
332 {
333 double noiselessActivation = computeNoiselessActivation (
334 preferredDirection,
335 actualDirection );
336
337 double probability = random.nextDouble ( );
338
339 double noisyActivation = Math.ceil (
340 normInv (
341 probability,
342 noiselessActivation,
343 RATE_DEVIATION ) );
344
345 if ( noisyActivation <= 0.0 )
346 {
347 noisyActivation = 0.0;
348 }
349
350 return noisyActivation;
351 }
352
353 public static double normInv (
354 double probability,
355 double mean,
356 double standardDeviation )
357 throws MathException
358 //////////////////////////////////////////////////////////////////////
359 {
360 NormalDistribution normalDistribution
361 = distributionFactory.createNormalDistribution (
362 mean, standardDeviation );
363
364 return
365 normalDistribution.inverseCumulativeProbability ( probability );
366 }
367
368 /*********************************************************************
369 * Predicts the direction for this trial.
370 * Knowing only the preferred directions for each neuron and their
371 * current activations, predict the direction for this trial.
372 * This method is the heart of the algorithm.
373 *********************************************************************/
374 public static double predictDirection (
375 double [ ] preferredDirections,
376 double [ ] noisyActivations )
377 //////////////////////////////////////////////////////////////////////
378 {
379 double sumCosComponents = 0.0;
380
381 double sumSinComponents = 0.0;
382
383 for ( int unit = 0; unit < NEURONS; unit++ )
384 {
385 // weight the activation of each neuron by their cartesian
386 // components along their preferred directions
387
388 double cosComponent
389 = Math.cos ( preferredDirections [ unit ] )
390 * noisyActivations [ unit ];
391
392 double sinComponent
393 = Math.sin ( preferredDirections [ unit ] )
394 * noisyActivations [ unit ];
395
396 // sum the weighted contributions from each neuron for this trial
397
398 sumCosComponents += cosComponent;
399
400 sumSinComponents += sinComponent;
401 }
402
403 // Here is the predicted direction.
404
405 double predictedDirection = Math.atan2 (
406 sumSinComponents,
407 sumCosComponents );
408
409 // return a non-negative angle for the answer
410
411 if ( predictedDirection < 0.0 )
412 {
413 predictedDirection += 2.0 * Math.PI;
414 }
415
416 // Note the arctan adjustment is not necessary since
417 // the atan2(y,x) function was used.
418
419 return predictedDirection;
420 }
421
422 //////////////////////////////////////////////////////////////////////
423 //////////////////////////////////////////////////////////////////////
424 }