Context Navigation

source: trunk/sources/ALGLIB/mlpbase.cs @ 2575

Visit:

Last change on this file since 2575 was 2563, checked in by gkronber, 15 years ago
Updated ALGLIB to latest version. #751 (Plugin for for data-modeling with ANN (integrated into CEDMA))
File size: 126.6 KB

Line
1	/*************************************************************************
2	Copyright (c) 2007-2008, Sergey Bochkanov (ALGLIB project).
3
4	>>> SOURCE LICENSE >>>
5	This program is free software; you can redistribute it and/or modify
6	it under the terms of the GNU General Public License as published by
7	the Free Software Foundation (www.fsf.org); either version 2 of the
8	License, or (at your option) any later version.
9
10	This program is distributed in the hope that it will be useful,
11	but WITHOUT ANY WARRANTY; without even the implied warranty of
12	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13	GNU General Public License for more details.
14
15	A copy of the GNU General Public License is available at
16	http://www.fsf.org/licensing/licenses
17
18	>>> END OF LICENSE >>>
19	*************************************************************************/
20
21	using System;
22
23	namespace alglib
24	{
25	public class mlpbase
26	{
27	public struct multilayerperceptron
28	{
29	public int[] structinfo;
30	public double[] weights;
31	public double[] columnmeans;
32	public double[] columnsigmas;
33	public double[] neurons;
34	public double[] dfdnet;
35	public double[] derror;
36	public double[] x;
37	public double[] y;
38	public double[,] chunks;
39	public double[] nwbuf;
40	};
41
42
43
44
45	public const int mlpvnum = 7;
46	public const int nfieldwidth = 4;
47	public const int chunksize = 32;
48
49
50	/*************************************************************************
51	Creates neural network with NIn inputs, NOut outputs, without hidden
52	layers, with linear output layer. Network weights are filled with small
53	random values.
54
55	-- ALGLIB --
56	Copyright 04.11.2007 by Bochkanov Sergey
57	*************************************************************************/
58	public static void mlpcreate0(int nin,
59	int nout,
60	ref multilayerperceptron network)
61	{
62	int[] lsizes = new int[0];
63	int[] ltypes = new int[0];
64	int[] lconnfirst = new int[0];
65	int[] lconnlast = new int[0];
66	int layerscount = 0;
67	int lastproc = 0;
68
69	layerscount = 1+2;
70
71	//
72	// Allocate arrays
73	//
74	lsizes = new int[layerscount-1+1];
75	ltypes = new int[layerscount-1+1];
76	lconnfirst = new int[layerscount-1+1];
77	lconnlast = new int[layerscount-1+1];
78
79	//
80	// Layers
81	//
82	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
83	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
84
85	//
86	// Create
87	//
88	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
89	}
90
91
92	/*************************************************************************
93	Same as MLPCreate0, but with one hidden layer (NHid neurons) with
94	non-linear activation function. Output layer is linear.
95
96	-- ALGLIB --
97	Copyright 04.11.2007 by Bochkanov Sergey
98	*************************************************************************/
99	public static void mlpcreate1(int nin,
100	int nhid,
101	int nout,
102	ref multilayerperceptron network)
103	{
104	int[] lsizes = new int[0];
105	int[] ltypes = new int[0];
106	int[] lconnfirst = new int[0];
107	int[] lconnlast = new int[0];
108	int layerscount = 0;
109	int lastproc = 0;
110
111	layerscount = 1+3+2;
112
113	//
114	// Allocate arrays
115	//
116	lsizes = new int[layerscount-1+1];
117	ltypes = new int[layerscount-1+1];
118	lconnfirst = new int[layerscount-1+1];
119	lconnlast = new int[layerscount-1+1];
120
121	//
122	// Layers
123	//
124	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
125	addbiasedsummatorlayer(nhid, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
126	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
127	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
128
129	//
130	// Create
131	//
132	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
133	}
134
135
136	/*************************************************************************
137	Same as MLPCreate0, but with two hidden layers (NHid1 and NHid2 neurons)
138	with non-linear activation function. Output layer is linear.
139	$ALL
140
141	-- ALGLIB --
142	Copyright 04.11.2007 by Bochkanov Sergey
143	*************************************************************************/
144	public static void mlpcreate2(int nin,
145	int nhid1,
146	int nhid2,
147	int nout,
148	ref multilayerperceptron network)
149	{
150	int[] lsizes = new int[0];
151	int[] ltypes = new int[0];
152	int[] lconnfirst = new int[0];
153	int[] lconnlast = new int[0];
154	int layerscount = 0;
155	int lastproc = 0;
156
157	layerscount = 1+3+3+2;
158
159	//
160	// Allocate arrays
161	//
162	lsizes = new int[layerscount-1+1];
163	ltypes = new int[layerscount-1+1];
164	lconnfirst = new int[layerscount-1+1];
165	lconnlast = new int[layerscount-1+1];
166
167	//
168	// Layers
169	//
170	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
171	addbiasedsummatorlayer(nhid1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
172	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
173	addbiasedsummatorlayer(nhid2, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
174	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
175	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
176
177	//
178	// Create
179	//
180	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
181	}
182
183
184	/*************************************************************************
185	Creates neural network with NIn inputs, NOut outputs, without hidden
186	layers with non-linear output layer. Network weights are filled with small
187	random values.
188
189	Activation function of the output layer takes values:
190
191	(B, +INF), if D>=0
192
193	or
194
195	(-INF, B), if D<0.
196
197
198	-- ALGLIB --
199	Copyright 30.03.2008 by Bochkanov Sergey
200	*************************************************************************/
201	public static void mlpcreateb0(int nin,
202	int nout,
203	double b,
204	double d,
205	ref multilayerperceptron network)
206	{
207	int[] lsizes = new int[0];
208	int[] ltypes = new int[0];
209	int[] lconnfirst = new int[0];
210	int[] lconnlast = new int[0];
211	int layerscount = 0;
212	int lastproc = 0;
213	int i = 0;
214
215	layerscount = 1+3;
216	if( (double)(d)>=(double)(0) )
217	{
218	d = 1;
219	}
220	else
221	{
222	d = -1;
223	}
224
225	//
226	// Allocate arrays
227	//
228	lsizes = new int[layerscount-1+1];
229	ltypes = new int[layerscount-1+1];
230	lconnfirst = new int[layerscount-1+1];
231	lconnlast = new int[layerscount-1+1];
232
233	//
234	// Layers
235	//
236	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
237	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
238	addactivationlayer(3, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
239
240	//
241	// Create
242	//
243	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
244
245	//
246	// Turn on ouputs shift/scaling.
247	//
248	for(i=nin; i<=nin+nout-1; i++)
249	{
250	network.columnmeans[i] = b;
251	network.columnsigmas[i] = d;
252	}
253	}
254
255
256	/*************************************************************************
257	Same as MLPCreateB0 but with non-linear hidden layer.
258
259	-- ALGLIB --
260	Copyright 30.03.2008 by Bochkanov Sergey
261	*************************************************************************/
262	public static void mlpcreateb1(int nin,
263	int nhid,
264	int nout,
265	double b,
266	double d,
267	ref multilayerperceptron network)
268	{
269	int[] lsizes = new int[0];
270	int[] ltypes = new int[0];
271	int[] lconnfirst = new int[0];
272	int[] lconnlast = new int[0];
273	int layerscount = 0;
274	int lastproc = 0;
275	int i = 0;
276
277	layerscount = 1+3+3;
278	if( (double)(d)>=(double)(0) )
279	{
280	d = 1;
281	}
282	else
283	{
284	d = -1;
285	}
286
287	//
288	// Allocate arrays
289	//
290	lsizes = new int[layerscount-1+1];
291	ltypes = new int[layerscount-1+1];
292	lconnfirst = new int[layerscount-1+1];
293	lconnlast = new int[layerscount-1+1];
294
295	//
296	// Layers
297	//
298	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
299	addbiasedsummatorlayer(nhid, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
300	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
301	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
302	addactivationlayer(3, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
303
304	//
305	// Create
306	//
307	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
308
309	//
310	// Turn on ouputs shift/scaling.
311	//
312	for(i=nin; i<=nin+nout-1; i++)
313	{
314	network.columnmeans[i] = b;
315	network.columnsigmas[i] = d;
316	}
317	}
318
319
320	/*************************************************************************
321	Same as MLPCreateB0 but with two non-linear hidden layers.
322
323	-- ALGLIB --
324	Copyright 30.03.2008 by Bochkanov Sergey
325	*************************************************************************/
326	public static void mlpcreateb2(int nin,
327	int nhid1,
328	int nhid2,
329	int nout,
330	double b,
331	double d,
332	ref multilayerperceptron network)
333	{
334	int[] lsizes = new int[0];
335	int[] ltypes = new int[0];
336	int[] lconnfirst = new int[0];
337	int[] lconnlast = new int[0];
338	int layerscount = 0;
339	int lastproc = 0;
340	int i = 0;
341
342	layerscount = 1+3+3+3;
343	if( (double)(d)>=(double)(0) )
344	{
345	d = 1;
346	}
347	else
348	{
349	d = -1;
350	}
351
352	//
353	// Allocate arrays
354	//
355	lsizes = new int[layerscount-1+1];
356	ltypes = new int[layerscount-1+1];
357	lconnfirst = new int[layerscount-1+1];
358	lconnlast = new int[layerscount-1+1];
359
360	//
361	// Layers
362	//
363	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
364	addbiasedsummatorlayer(nhid1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
365	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
366	addbiasedsummatorlayer(nhid2, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
367	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
368	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
369	addactivationlayer(3, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
370
371	//
372	// Create
373	//
374	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
375
376	//
377	// Turn on ouputs shift/scaling.
378	//
379	for(i=nin; i<=nin+nout-1; i++)
380	{
381	network.columnmeans[i] = b;
382	network.columnsigmas[i] = d;
383	}
384	}
385
386
387	/*************************************************************************
388	Creates neural network with NIn inputs, NOut outputs, without hidden
389	layers with non-linear output layer. Network weights are filled with small
390	random values. Activation function of the output layer takes values [A,B].
391
392	-- ALGLIB --
393	Copyright 30.03.2008 by Bochkanov Sergey
394	*************************************************************************/
395	public static void mlpcreater0(int nin,
396	int nout,
397	double a,
398	double b,
399	ref multilayerperceptron network)
400	{
401	int[] lsizes = new int[0];
402	int[] ltypes = new int[0];
403	int[] lconnfirst = new int[0];
404	int[] lconnlast = new int[0];
405	int layerscount = 0;
406	int lastproc = 0;
407	int i = 0;
408
409	layerscount = 1+3;
410
411	//
412	// Allocate arrays
413	//
414	lsizes = new int[layerscount-1+1];
415	ltypes = new int[layerscount-1+1];
416	lconnfirst = new int[layerscount-1+1];
417	lconnlast = new int[layerscount-1+1];
418
419	//
420	// Layers
421	//
422	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
423	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
424	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
425
426	//
427	// Create
428	//
429	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
430
431	//
432	// Turn on outputs shift/scaling.
433	//
434	for(i=nin; i<=nin+nout-1; i++)
435	{
436	network.columnmeans[i] = 0.5*(a+b);
437	network.columnsigmas[i] = 0.5*(a-b);
438	}
439	}
440
441
442	/*************************************************************************
443	Same as MLPCreateR0, but with non-linear hidden layer.
444
445	-- ALGLIB --
446	Copyright 30.03.2008 by Bochkanov Sergey
447	*************************************************************************/
448	public static void mlpcreater1(int nin,
449	int nhid,
450	int nout,
451	double a,
452	double b,
453	ref multilayerperceptron network)
454	{
455	int[] lsizes = new int[0];
456	int[] ltypes = new int[0];
457	int[] lconnfirst = new int[0];
458	int[] lconnlast = new int[0];
459	int layerscount = 0;
460	int lastproc = 0;
461	int i = 0;
462
463	layerscount = 1+3+3;
464
465	//
466	// Allocate arrays
467	//
468	lsizes = new int[layerscount-1+1];
469	ltypes = new int[layerscount-1+1];
470	lconnfirst = new int[layerscount-1+1];
471	lconnlast = new int[layerscount-1+1];
472
473	//
474	// Layers
475	//
476	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
477	addbiasedsummatorlayer(nhid, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
478	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
479	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
480	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
481
482	//
483	// Create
484	//
485	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
486
487	//
488	// Turn on outputs shift/scaling.
489	//
490	for(i=nin; i<=nin+nout-1; i++)
491	{
492	network.columnmeans[i] = 0.5*(a+b);
493	network.columnsigmas[i] = 0.5*(a-b);
494	}
495	}
496
497
498	/*************************************************************************
499	Same as MLPCreateR0, but with two non-linear hidden layers.
500
501	-- ALGLIB --
502	Copyright 30.03.2008 by Bochkanov Sergey
503	*************************************************************************/
504	public static void mlpcreater2(int nin,
505	int nhid1,
506	int nhid2,
507	int nout,
508	double a,
509	double b,
510	ref multilayerperceptron network)
511	{
512	int[] lsizes = new int[0];
513	int[] ltypes = new int[0];
514	int[] lconnfirst = new int[0];
515	int[] lconnlast = new int[0];
516	int layerscount = 0;
517	int lastproc = 0;
518	int i = 0;
519
520	layerscount = 1+3+3+3;
521
522	//
523	// Allocate arrays
524	//
525	lsizes = new int[layerscount-1+1];
526	ltypes = new int[layerscount-1+1];
527	lconnfirst = new int[layerscount-1+1];
528	lconnlast = new int[layerscount-1+1];
529
530	//
531	// Layers
532	//
533	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
534	addbiasedsummatorlayer(nhid1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
535	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
536	addbiasedsummatorlayer(nhid2, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
537	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
538	addbiasedsummatorlayer(nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
539	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
540
541	//
542	// Create
543	//
544	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, false, ref network);
545
546	//
547	// Turn on outputs shift/scaling.
548	//
549	for(i=nin; i<=nin+nout-1; i++)
550	{
551	network.columnmeans[i] = 0.5*(a+b);
552	network.columnsigmas[i] = 0.5*(a-b);
553	}
554	}
555
556
557	/*************************************************************************
558	Creates classifier network with NIn inputs and NOut possible classes.
559	Network contains no hidden layers and linear output layer with SOFTMAX-
560	normalization (so outputs sums up to 1.0 and converge to posterior
561	probabilities).
562
563	-- ALGLIB --
564	Copyright 04.11.2007 by Bochkanov Sergey
565	*************************************************************************/
566	public static void mlpcreatec0(int nin,
567	int nout,
568	ref multilayerperceptron network)
569	{
570	int[] lsizes = new int[0];
571	int[] ltypes = new int[0];
572	int[] lconnfirst = new int[0];
573	int[] lconnlast = new int[0];
574	int layerscount = 0;
575	int lastproc = 0;
576
577	System.Diagnostics.Debug.Assert(nout>=2, "MLPCreateC0: NOut<2!");
578	layerscount = 1+2+1;
579
580	//
581	// Allocate arrays
582	//
583	lsizes = new int[layerscount-1+1];
584	ltypes = new int[layerscount-1+1];
585	lconnfirst = new int[layerscount-1+1];
586	lconnlast = new int[layerscount-1+1];
587
588	//
589	// Layers
590	//
591	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
592	addbiasedsummatorlayer(nout-1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
593	addzerolayer(ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
594
595	//
596	// Create
597	//
598	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, true, ref network);
599	}
600
601
602	/*************************************************************************
603	Same as MLPCreateC0, but with one non-linear hidden layer.
604
605	-- ALGLIB --
606	Copyright 04.11.2007 by Bochkanov Sergey
607	*************************************************************************/
608	public static void mlpcreatec1(int nin,
609	int nhid,
610	int nout,
611	ref multilayerperceptron network)
612	{
613	int[] lsizes = new int[0];
614	int[] ltypes = new int[0];
615	int[] lconnfirst = new int[0];
616	int[] lconnlast = new int[0];
617	int layerscount = 0;
618	int lastproc = 0;
619
620	System.Diagnostics.Debug.Assert(nout>=2, "MLPCreateC1: NOut<2!");
621	layerscount = 1+3+2+1;
622
623	//
624	// Allocate arrays
625	//
626	lsizes = new int[layerscount-1+1];
627	ltypes = new int[layerscount-1+1];
628	lconnfirst = new int[layerscount-1+1];
629	lconnlast = new int[layerscount-1+1];
630
631	//
632	// Layers
633	//
634	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
635	addbiasedsummatorlayer(nhid, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
636	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
637	addbiasedsummatorlayer(nout-1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
638	addzerolayer(ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
639
640	//
641	// Create
642	//
643	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, true, ref network);
644	}
645
646
647	/*************************************************************************
648	Same as MLPCreateC0, but with two non-linear hidden layers.
649
650	-- ALGLIB --
651	Copyright 04.11.2007 by Bochkanov Sergey
652	*************************************************************************/
653	public static void mlpcreatec2(int nin,
654	int nhid1,
655	int nhid2,
656	int nout,
657	ref multilayerperceptron network)
658	{
659	int[] lsizes = new int[0];
660	int[] ltypes = new int[0];
661	int[] lconnfirst = new int[0];
662	int[] lconnlast = new int[0];
663	int layerscount = 0;
664	int lastproc = 0;
665
666	System.Diagnostics.Debug.Assert(nout>=2, "MLPCreateC2: NOut<2!");
667	layerscount = 1+3+3+2+1;
668
669	//
670	// Allocate arrays
671	//
672	lsizes = new int[layerscount-1+1];
673	ltypes = new int[layerscount-1+1];
674	lconnfirst = new int[layerscount-1+1];
675	lconnlast = new int[layerscount-1+1];
676
677	//
678	// Layers
679	//
680	addinputlayer(nin, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
681	addbiasedsummatorlayer(nhid1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
682	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
683	addbiasedsummatorlayer(nhid2, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
684	addactivationlayer(1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
685	addbiasedsummatorlayer(nout-1, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
686	addzerolayer(ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, ref lastproc);
687
688	//
689	// Create
690	//
691	mlpcreate(nin, nout, ref lsizes, ref ltypes, ref lconnfirst, ref lconnlast, layerscount, true, ref network);
692	}
693
694
695	/*************************************************************************
696	Copying of neural network
697
698	INPUT PARAMETERS:
699	Network1 - original
700
701	OUTPUT PARAMETERS:
702	Network2 - copy
703
704	-- ALGLIB --
705	Copyright 04.11.2007 by Bochkanov Sergey
706	*************************************************************************/
707	public static void mlpcopy(ref multilayerperceptron network1,
708	ref multilayerperceptron network2)
709	{
710	int i = 0;
711	int ssize = 0;
712	int ntotal = 0;
713	int nin = 0;
714	int nout = 0;
715	int wcount = 0;
716	int i_ = 0;
717
718
719	//
720	// Unload info
721	//
722	ssize = network1.structinfo[0];
723	nin = network1.structinfo[1];
724	nout = network1.structinfo[2];
725	ntotal = network1.structinfo[3];
726	wcount = network1.structinfo[4];
727
728	//
729	// Allocate space
730	//
731	network2.structinfo = new int[ssize-1+1];
732	network2.weights = new double[wcount-1+1];
733	if( mlpissoftmax(ref network1) )
734	{
735	network2.columnmeans = new double[nin-1+1];
736	network2.columnsigmas = new double[nin-1+1];
737	}
738	else
739	{
740	network2.columnmeans = new double[nin+nout-1+1];
741	network2.columnsigmas = new double[nin+nout-1+1];
742	}
743	network2.neurons = new double[ntotal-1+1];
744	network2.chunks = new double[3*ntotal+1, chunksize-1+1];
745	network2.nwbuf = new double[Math.Max(wcount, 2*nout)-1+1];
746	network2.dfdnet = new double[ntotal-1+1];
747	network2.x = new double[nin-1+1];
748	network2.y = new double[nout-1+1];
749	network2.derror = new double[ntotal-1+1];
750
751	//
752	// Copy
753	//
754	for(i=0; i<=ssize-1; i++)
755	{
756	network2.structinfo[i] = network1.structinfo[i];
757	}
758	for(i_=0; i_<=wcount-1;i_++)
759	{
760	network2.weights[i_] = network1.weights[i_];
761	}
762	if( mlpissoftmax(ref network1) )
763	{
764	for(i_=0; i_<=nin-1;i_++)
765	{
766	network2.columnmeans[i_] = network1.columnmeans[i_];
767	}
768	for(i_=0; i_<=nin-1;i_++)
769	{
770	network2.columnsigmas[i_] = network1.columnsigmas[i_];
771	}
772	}
773	else
774	{
775	for(i_=0; i_<=nin+nout-1;i_++)
776	{
777	network2.columnmeans[i_] = network1.columnmeans[i_];
778	}
779	for(i_=0; i_<=nin+nout-1;i_++)
780	{
781	network2.columnsigmas[i_] = network1.columnsigmas[i_];
782	}
783	}
784	for(i_=0; i_<=ntotal-1;i_++)
785	{
786	network2.neurons[i_] = network1.neurons[i_];
787	}
788	for(i_=0; i_<=ntotal-1;i_++)
789	{
790	network2.dfdnet[i_] = network1.dfdnet[i_];
791	}
792	for(i_=0; i_<=nin-1;i_++)
793	{
794	network2.x[i_] = network1.x[i_];
795	}
796	for(i_=0; i_<=nout-1;i_++)
797	{
798	network2.y[i_] = network1.y[i_];
799	}
800	for(i_=0; i_<=ntotal-1;i_++)
801	{
802	network2.derror[i_] = network1.derror[i_];
803	}
804	}
805
806
807	/*************************************************************************
808	Serialization of MultiLayerPerceptron strucure
809
810	INPUT PARAMETERS:
811	Network - original
812
813	OUTPUT PARAMETERS:
814	RA - array of real numbers which stores network,
815	array[0..RLen-1]
816	RLen - RA lenght
817
818	-- ALGLIB --
819	Copyright 29.03.2008 by Bochkanov Sergey
820	*************************************************************************/
821	public static void mlpserialize(ref multilayerperceptron network,
822	ref double[] ra,
823	ref int rlen)
824	{
825	int i = 0;
826	int ssize = 0;
827	int ntotal = 0;
828	int nin = 0;
829	int nout = 0;
830	int wcount = 0;
831	int sigmalen = 0;
832	int offs = 0;
833	int i_ = 0;
834	int i1_ = 0;
835
836
837	//
838	// Unload info
839	//
840	ssize = network.structinfo[0];
841	nin = network.structinfo[1];
842	nout = network.structinfo[2];
843	ntotal = network.structinfo[3];
844	wcount = network.structinfo[4];
845	if( mlpissoftmax(ref network) )
846	{
847	sigmalen = nin;
848	}
849	else
850	{
851	sigmalen = nin+nout;
852	}
853
854	//
855	// RA format:
856	// LEN DESRC.
857	// 1 RLen
858	// 1 version (MLPVNum)
859	// 1 StructInfo size
860	// SSize StructInfo
861	// WCount Weights
862	// SigmaLen ColumnMeans
863	// SigmaLen ColumnSigmas
864	//
865	rlen = 3+ssize+wcount+2*sigmalen;
866	ra = new double[rlen-1+1];
867	ra[0] = rlen;
868	ra[1] = mlpvnum;
869	ra[2] = ssize;
870	offs = 3;
871	for(i=0; i<=ssize-1; i++)
872	{
873	ra[offs+i] = network.structinfo[i];
874	}
875	offs = offs+ssize;
876	i1_ = (0) - (offs);
877	for(i_=offs; i_<=offs+wcount-1;i_++)
878	{
879	ra[i_] = network.weights[i_+i1_];
880	}
881	offs = offs+wcount;
882	i1_ = (0) - (offs);
883	for(i_=offs; i_<=offs+sigmalen-1;i_++)
884	{
885	ra[i_] = network.columnmeans[i_+i1_];
886	}
887	offs = offs+sigmalen;
888	i1_ = (0) - (offs);
889	for(i_=offs; i_<=offs+sigmalen-1;i_++)
890	{
891	ra[i_] = network.columnsigmas[i_+i1_];
892	}
893	offs = offs+sigmalen;
894	}
895
896
897	/*************************************************************************
898	Unserialization of MultiLayerPerceptron strucure
899
900	INPUT PARAMETERS:
901	RA - real array which stores network
902
903	OUTPUT PARAMETERS:
904	Network - restored network
905
906	-- ALGLIB --
907	Copyright 29.03.2008 by Bochkanov Sergey
908	*************************************************************************/
909	public static void mlpunserialize(ref double[] ra,
910	ref multilayerperceptron network)
911	{
912	int i = 0;
913	int ssize = 0;
914	int ntotal = 0;
915	int nin = 0;
916	int nout = 0;
917	int wcount = 0;
918	int sigmalen = 0;
919	int offs = 0;
920	int i_ = 0;
921	int i1_ = 0;
922
923	System.Diagnostics.Debug.Assert((int)Math.Round(ra[1])==mlpvnum, "MLPUnserialize: incorrect array!");
924
925	//
926	// Unload StructInfo from IA
927	//
928	offs = 3;
929	ssize = (int)Math.Round(ra[2]);
930	network.structinfo = new int[ssize-1+1];
931	for(i=0; i<=ssize-1; i++)
932	{
933	network.structinfo[i] = (int)Math.Round(ra[offs+i]);
934	}
935	offs = offs+ssize;
936
937	//
938	// Unload info from StructInfo
939	//
940	ssize = network.structinfo[0];
941	nin = network.structinfo[1];
942	nout = network.structinfo[2];
943	ntotal = network.structinfo[3];
944	wcount = network.structinfo[4];
945	if( network.structinfo[6]==0 )
946	{
947	sigmalen = nin+nout;
948	}
949	else
950	{
951	sigmalen = nin;
952	}
953
954	//
955	// Allocate space for other fields
956	//
957	network.weights = new double[wcount-1+1];
958	network.columnmeans = new double[sigmalen-1+1];
959	network.columnsigmas = new double[sigmalen-1+1];
960	network.neurons = new double[ntotal-1+1];
961	network.chunks = new double[3*ntotal+1, chunksize-1+1];
962	network.nwbuf = new double[Math.Max(wcount, 2*nout)-1+1];
963	network.dfdnet = new double[ntotal-1+1];
964	network.x = new double[nin-1+1];
965	network.y = new double[nout-1+1];
966	network.derror = new double[ntotal-1+1];
967
968	//
969	// Copy parameters from RA
970	//
971	i1_ = (offs) - (0);
972	for(i_=0; i_<=wcount-1;i_++)
973	{
974	network.weights[i_] = ra[i_+i1_];
975	}
976	offs = offs+wcount;
977	i1_ = (offs) - (0);
978	for(i_=0; i_<=sigmalen-1;i_++)
979	{
980	network.columnmeans[i_] = ra[i_+i1_];
981	}
982	offs = offs+sigmalen;
983	i1_ = (offs) - (0);
984	for(i_=0; i_<=sigmalen-1;i_++)
985	{
986	network.columnsigmas[i_] = ra[i_+i1_];
987	}
988	offs = offs+sigmalen;
989	}
990
991
992	/*************************************************************************
993	Randomization of neural network weights
994
995	-- ALGLIB --
996	Copyright 06.11.2007 by Bochkanov Sergey
997	*************************************************************************/
998	public static void mlprandomize(ref multilayerperceptron network)
999	{
1000	int i = 0;
1001	int nin = 0;
1002	int nout = 0;
1003	int wcount = 0;
1004
1005	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1006	for(i=0; i<=wcount-1; i++)
1007	{
1008	network.weights[i] = AP.Math.RandomReal()-0.5;
1009	}
1010	}
1011
1012
1013	/*************************************************************************
1014	Randomization of neural network weights and standartisator
1015
1016	-- ALGLIB --
1017	Copyright 10.03.2008 by Bochkanov Sergey
1018	*************************************************************************/
1019	public static void mlprandomizefull(ref multilayerperceptron network)
1020	{
1021	int i = 0;
1022	int nin = 0;
1023	int nout = 0;
1024	int wcount = 0;
1025	int ntotal = 0;
1026	int istart = 0;
1027	int offs = 0;
1028	int ntype = 0;
1029
1030	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1031	ntotal = network.structinfo[3];
1032	istart = network.structinfo[5];
1033
1034	//
1035	// Process network
1036	//
1037	for(i=0; i<=wcount-1; i++)
1038	{
1039	network.weights[i] = AP.Math.RandomReal()-0.5;
1040	}
1041	for(i=0; i<=nin-1; i++)
1042	{
1043	network.columnmeans[i] = 2*AP.Math.RandomReal()-1;
1044	network.columnsigmas[i] = 1.5*AP.Math.RandomReal()+0.5;
1045	}
1046	if( !mlpissoftmax(ref network) )
1047	{
1048	for(i=0; i<=nout-1; i++)
1049	{
1050	offs = istart+(ntotal-nout+i)*nfieldwidth;
1051	ntype = network.structinfo[offs+0];
1052	if( ntype==0 )
1053	{
1054
1055	//
1056	// Shifts are changed only for linear outputs neurons
1057	//
1058	network.columnmeans[nin+i] = 2*AP.Math.RandomReal()-1;
1059	}
1060	if( ntype==0 \| ntype==3 )
1061	{
1062
1063	//
1064	// Scales are changed only for linear or bounded outputs neurons.
1065	// Note that scale randomization preserves sign.
1066	//
1067	network.columnsigmas[nin+i] = Math.Sign(network.columnsigmas[nin+i])(1.5AP.Math.RandomReal()+0.5);
1068	}
1069	}
1070	}
1071	}
1072
1073
1074	/*************************************************************************
1075	Internal subroutine.
1076
1077	-- ALGLIB --
1078	Copyright 30.03.2008 by Bochkanov Sergey
1079	*************************************************************************/
1080	public static void mlpinitpreprocessor(ref multilayerperceptron network,
1081	ref double[,] xy,
1082	int ssize)
1083	{
1084	int i = 0;
1085	int j = 0;
1086	int jmax = 0;
1087	int nin = 0;
1088	int nout = 0;
1089	int wcount = 0;
1090	int ntotal = 0;
1091	int istart = 0;
1092	int offs = 0;
1093	int ntype = 0;
1094	double[] means = new double[0];
1095	double[] sigmas = new double[0];
1096	double s = 0;
1097
1098	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1099	ntotal = network.structinfo[3];
1100	istart = network.structinfo[5];
1101
1102	//
1103	// Means/Sigmas
1104	//
1105	if( mlpissoftmax(ref network) )
1106	{
1107	jmax = nin-1;
1108	}
1109	else
1110	{
1111	jmax = nin+nout-1;
1112	}
1113	means = new double[jmax+1];
1114	sigmas = new double[jmax+1];
1115	for(j=0; j<=jmax; j++)
1116	{
1117	means[j] = 0;
1118	for(i=0; i<=ssize-1; i++)
1119	{
1120	means[j] = means[j]+xy[i,j];
1121	}
1122	means[j] = means[j]/ssize;
1123	sigmas[j] = 0;
1124	for(i=0; i<=ssize-1; i++)
1125	{
1126	sigmas[j] = sigmas[j]+AP.Math.Sqr(xy[i,j]-means[j]);
1127	}
1128	sigmas[j] = Math.Sqrt(sigmas[j]/ssize);
1129	}
1130
1131	//
1132	// Inputs
1133	//
1134	for(i=0; i<=nin-1; i++)
1135	{
1136	network.columnmeans[i] = means[i];
1137	network.columnsigmas[i] = sigmas[i];
1138	if( (double)(network.columnsigmas[i])==(double)(0) )
1139	{
1140	network.columnsigmas[i] = 1;
1141	}
1142	}
1143
1144	//
1145	// Outputs
1146	//
1147	if( !mlpissoftmax(ref network) )
1148	{
1149	for(i=0; i<=nout-1; i++)
1150	{
1151	offs = istart+(ntotal-nout+i)*nfieldwidth;
1152	ntype = network.structinfo[offs+0];
1153
1154	//
1155	// Linear outputs
1156	//
1157	if( ntype==0 )
1158	{
1159	network.columnmeans[nin+i] = means[nin+i];
1160	network.columnsigmas[nin+i] = sigmas[nin+i];
1161	if( (double)(network.columnsigmas[nin+i])==(double)(0) )
1162	{
1163	network.columnsigmas[nin+i] = 1;
1164	}
1165	}
1166
1167	//
1168	// Bounded outputs (half-interval)
1169	//
1170	if( ntype==3 )
1171	{
1172	s = means[nin+i]-network.columnmeans[nin+i];
1173	if( (double)(s)==(double)(0) )
1174	{
1175	s = Math.Sign(network.columnsigmas[nin+i]);
1176	}
1177	if( (double)(s)==(double)(0) )
1178	{
1179	s = 1.0;
1180	}
1181	network.columnsigmas[nin+i] = Math.Sign(network.columnsigmas[nin+i])*Math.Abs(s);
1182	if( (double)(network.columnsigmas[nin+i])==(double)(0) )
1183	{
1184	network.columnsigmas[nin+i] = 1;
1185	}
1186	}
1187	}
1188	}
1189	}
1190
1191
1192	/*************************************************************************
1193	Returns information about initialized network: number of inputs, outputs,
1194	weights.
1195
1196	-- ALGLIB --
1197	Copyright 04.11.2007 by Bochkanov Sergey
1198	*************************************************************************/
1199	public static void mlpproperties(ref multilayerperceptron network,
1200	ref int nin,
1201	ref int nout,
1202	ref int wcount)
1203	{
1204	nin = network.structinfo[1];
1205	nout = network.structinfo[2];
1206	wcount = network.structinfo[4];
1207	}
1208
1209
1210	/*************************************************************************
1211	Tells whether network is SOFTMAX-normalized (i.e. classifier) or not.
1212
1213	-- ALGLIB --
1214	Copyright 04.11.2007 by Bochkanov Sergey
1215	*************************************************************************/
1216	public static bool mlpissoftmax(ref multilayerperceptron network)
1217	{
1218	bool result = new bool();
1219
1220	result = network.structinfo[6]==1;
1221	return result;
1222	}
1223
1224
1225	/*************************************************************************
1226	Procesing
1227
1228	INPUT PARAMETERS:
1229	Network - neural network
1230	X - input vector, array[0..NIn-1].
1231
1232	OUTPUT PARAMETERS:
1233	Y - result. Regression estimate when solving regression task,
1234	vector of posterior probabilities for classification task.
1235	Subroutine does not allocate memory for this vector, it is
1236	responsibility of a caller to allocate it. Array must be
1237	at least [0..NOut-1].
1238
1239	-- ALGLIB --
1240	Copyright 04.11.2007 by Bochkanov Sergey
1241	*************************************************************************/
1242	public static void mlpprocess(ref multilayerperceptron network,
1243	ref double[] x,
1244	ref double[] y)
1245	{
1246	mlpinternalprocessvector(ref network.structinfo, ref network.weights, ref network.columnmeans, ref network.columnsigmas, ref network.neurons, ref network.dfdnet, ref x, ref y);
1247	}
1248
1249
1250	/*************************************************************************
1251	Error function for neural network, internal subroutine.
1252
1253	-- ALGLIB --
1254	Copyright 04.11.2007 by Bochkanov Sergey
1255	*************************************************************************/
1256	public static double mlperror(ref multilayerperceptron network,
1257	ref double[,] xy,
1258	int ssize)
1259	{
1260	double result = 0;
1261	int i = 0;
1262	int k = 0;
1263	int nin = 0;
1264	int nout = 0;
1265	int wcount = 0;
1266	double e = 0;
1267	int i_ = 0;
1268	int i1_ = 0;
1269
1270	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1271	result = 0;
1272	for(i=0; i<=ssize-1; i++)
1273	{
1274	for(i_=0; i_<=nin-1;i_++)
1275	{
1276	network.x[i_] = xy[i,i_];
1277	}
1278	mlpprocess(ref network, ref network.x, ref network.y);
1279	if( mlpissoftmax(ref network) )
1280	{
1281
1282	//
1283	// class labels outputs
1284	//
1285	k = (int)Math.Round(xy[i,nin]);
1286	if( k>=0 & k<nout )
1287	{
1288	network.y[k] = network.y[k]-1;
1289	}
1290	}
1291	else
1292	{
1293
1294	//
1295	// real outputs
1296	//
1297	i1_ = (nin) - (0);
1298	for(i_=0; i_<=nout-1;i_++)
1299	{
1300	network.y[i_] = network.y[i_] - xy[i,i_+i1_];
1301	}
1302	}
1303	e = 0.0;
1304	for(i_=0; i_<=nout-1;i_++)
1305	{
1306	e += network.y[i_]*network.y[i_];
1307	}
1308	result = result+e/2;
1309	}
1310	return result;
1311	}
1312
1313
1314	/*************************************************************************
1315	Natural error function for neural network, internal subroutine.
1316
1317	-- ALGLIB --
1318	Copyright 04.11.2007 by Bochkanov Sergey
1319	*************************************************************************/
1320	public static double mlperrorn(ref multilayerperceptron network,
1321	ref double[,] xy,
1322	int ssize)
1323	{
1324	double result = 0;
1325	int i = 0;
1326	int j = 0;
1327	int k = 0;
1328	int nin = 0;
1329	int nout = 0;
1330	int wcount = 0;
1331	double e = 0;
1332	double t = 0;
1333	double z = 0;
1334	int i_ = 0;
1335	int i1_ = 0;
1336
1337	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1338	result = 0;
1339	for(i=0; i<=ssize-1; i++)
1340	{
1341
1342	//
1343	// Process vector
1344	//
1345	for(i_=0; i_<=nin-1;i_++)
1346	{
1347	network.x[i_] = xy[i,i_];
1348	}
1349	mlpprocess(ref network, ref network.x, ref network.y);
1350
1351	//
1352	// Update error function
1353	//
1354	if( network.structinfo[6]==0 )
1355	{
1356
1357	//
1358	// Least squares error function
1359	//
1360	i1_ = (nin) - (0);
1361	for(i_=0; i_<=nout-1;i_++)
1362	{
1363	network.y[i_] = network.y[i_] - xy[i,i_+i1_];
1364	}
1365	e = 0.0;
1366	for(i_=0; i_<=nout-1;i_++)
1367	{
1368	e += network.y[i_]*network.y[i_];
1369	}
1370	result = result+e/2;
1371	}
1372	else
1373	{
1374
1375	//
1376	// Cross-entropy error function
1377	//
1378	k = (int)Math.Round(xy[i,nin]);
1379	if( k>=0 & k<nout )
1380	{
1381	result = result+safecrossentropy(1, network.y[k]);
1382	}
1383	}
1384	}
1385	return result;
1386	}
1387
1388
1389	/*************************************************************************
1390	Classification error
1391
1392	-- ALGLIB --
1393	Copyright 04.11.2007 by Bochkanov Sergey
1394	*************************************************************************/
1395	public static int mlpclserror(ref multilayerperceptron network,
1396	ref double[,] xy,
1397	int ssize)
1398	{
1399	int result = 0;
1400	int i = 0;
1401	int j = 0;
1402	int nin = 0;
1403	int nout = 0;
1404	int wcount = 0;
1405	double[] workx = new double[0];
1406	double[] worky = new double[0];
1407	int nn = 0;
1408	int ns = 0;
1409	int nmax = 0;
1410	int i_ = 0;
1411
1412	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1413	workx = new double[nin-1+1];
1414	worky = new double[nout-1+1];
1415	result = 0;
1416	for(i=0; i<=ssize-1; i++)
1417	{
1418
1419	//
1420	// Process
1421	//
1422	for(i_=0; i_<=nin-1;i_++)
1423	{
1424	workx[i_] = xy[i,i_];
1425	}
1426	mlpprocess(ref network, ref workx, ref worky);
1427
1428	//
1429	// Network version of the answer
1430	//
1431	nmax = 0;
1432	for(j=0; j<=nout-1; j++)
1433	{
1434	if( (double)(worky[j])>(double)(worky[nmax]) )
1435	{
1436	nmax = j;
1437	}
1438	}
1439	nn = nmax;
1440
1441	//
1442	// Right answer
1443	//
1444	if( mlpissoftmax(ref network) )
1445	{
1446	ns = (int)Math.Round(xy[i,nin]);
1447	}
1448	else
1449	{
1450	nmax = 0;
1451	for(j=0; j<=nout-1; j++)
1452	{
1453	if( (double)(xy[i,nin+j])>(double)(xy[i,nin+nmax]) )
1454	{
1455	nmax = j;
1456	}
1457	}
1458	ns = nmax;
1459	}
1460
1461	//
1462	// compare
1463	//
1464	if( nn!=ns )
1465	{
1466	result = result+1;
1467	}
1468	}
1469	return result;
1470	}
1471
1472
1473	/*************************************************************************
1474	Relative classification error on the test set
1475
1476	INPUT PARAMETERS:
1477	Network - network
1478	XY - test set
1479	NPoints - test set size
1480
1481	RESULT:
1482	percent of incorrectly classified cases. Works both for
1483	classifier networks and general purpose networks used as
1484	classifiers.
1485
1486	-- ALGLIB --
1487	Copyright 25.12.2008 by Bochkanov Sergey
1488	*************************************************************************/
1489	public static double mlprelclserror(ref multilayerperceptron network,
1490	ref double[,] xy,
1491	int npoints)
1492	{
1493	double result = 0;
1494
1495	result = (double)(mlpclserror(ref network, ref xy, npoints))/(double)(npoints);
1496	return result;
1497	}
1498
1499
1500	/*************************************************************************
1501	Average cross-entropy (in bits per element) on the test set
1502
1503	INPUT PARAMETERS:
1504	Network - neural network
1505	XY - test set
1506	NPoints - test set size
1507
1508	RESULT:
1509	CrossEntropy/(NPoints*LN(2)).
1510	Zero if network solves regression task.
1511
1512	-- ALGLIB --
1513	Copyright 08.01.2009 by Bochkanov Sergey
1514	*************************************************************************/
1515	public static double mlpavgce(ref multilayerperceptron network,
1516	ref double[,] xy,
1517	int npoints)
1518	{
1519	double result = 0;
1520	int nin = 0;
1521	int nout = 0;
1522	int wcount = 0;
1523
1524	if( mlpissoftmax(ref network) )
1525	{
1526	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1527	result = mlperrorn(ref network, ref xy, npoints)/(npoints*Math.Log(2));
1528	}
1529	else
1530	{
1531	result = 0;
1532	}
1533	return result;
1534	}
1535
1536
1537	/*************************************************************************
1538	RMS error on the test set
1539
1540	INPUT PARAMETERS:
1541	Network - neural network
1542	XY - test set
1543	NPoints - test set size
1544
1545	RESULT:
1546	root mean square error.
1547	Its meaning for regression task is obvious. As for
1548	classification task, RMS error means error when estimating posterior
1549	probabilities.
1550
1551	-- ALGLIB --
1552	Copyright 04.11.2007 by Bochkanov Sergey
1553	*************************************************************************/
1554	public static double mlprmserror(ref multilayerperceptron network,
1555	ref double[,] xy,
1556	int npoints)
1557	{
1558	double result = 0;
1559	int nin = 0;
1560	int nout = 0;
1561	int wcount = 0;
1562
1563	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1564	result = Math.Sqrt(2mlperror(ref network, ref xy, npoints)/(npointsnout));
1565	return result;
1566	}
1567
1568
1569	/*************************************************************************
1570	Average error on the test set
1571
1572	INPUT PARAMETERS:
1573	Network - neural network
1574	XY - test set
1575	NPoints - test set size
1576
1577	RESULT:
1578	Its meaning for regression task is obvious. As for
1579	classification task, it means average error when estimating posterior
1580	probabilities.
1581
1582	-- ALGLIB --
1583	Copyright 11.03.2008 by Bochkanov Sergey
1584	*************************************************************************/
1585	public static double mlpavgerror(ref multilayerperceptron network,
1586	ref double[,] xy,
1587	int npoints)
1588	{
1589	double result = 0;
1590	int i = 0;
1591	int j = 0;
1592	int k = 0;
1593	int nin = 0;
1594	int nout = 0;
1595	int wcount = 0;
1596	int i_ = 0;
1597
1598	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1599	result = 0;
1600	for(i=0; i<=npoints-1; i++)
1601	{
1602	for(i_=0; i_<=nin-1;i_++)
1603	{
1604	network.x[i_] = xy[i,i_];
1605	}
1606	mlpprocess(ref network, ref network.x, ref network.y);
1607	if( mlpissoftmax(ref network) )
1608	{
1609
1610	//
1611	// class labels
1612	//
1613	k = (int)Math.Round(xy[i,nin]);
1614	for(j=0; j<=nout-1; j++)
1615	{
1616	if( j==k )
1617	{
1618	result = result+Math.Abs(1-network.y[j]);
1619	}
1620	else
1621	{
1622	result = result+Math.Abs(network.y[j]);
1623	}
1624	}
1625	}
1626	else
1627	{
1628
1629	//
1630	// real outputs
1631	//
1632	for(j=0; j<=nout-1; j++)
1633	{
1634	result = result+Math.Abs(xy[i,nin+j]-network.y[j]);
1635	}
1636	}
1637	}
1638	result = result/(npoints*nout);
1639	return result;
1640	}
1641
1642
1643	/*************************************************************************
1644	Average relative error on the test set
1645
1646	INPUT PARAMETERS:
1647	Network - neural network
1648	XY - test set
1649	NPoints - test set size
1650
1651	RESULT:
1652	Its meaning for regression task is obvious. As for
1653	classification task, it means average relative error when estimating
1654	posterior probability of belonging to the correct class.
1655
1656	-- ALGLIB --
1657	Copyright 11.03.2008 by Bochkanov Sergey
1658	*************************************************************************/
1659	public static double mlpavgrelerror(ref multilayerperceptron network,
1660	ref double[,] xy,
1661	int npoints)
1662	{
1663	double result = 0;
1664	int i = 0;
1665	int j = 0;
1666	int k = 0;
1667	int lk = 0;
1668	int nin = 0;
1669	int nout = 0;
1670	int wcount = 0;
1671	int i_ = 0;
1672
1673	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1674	result = 0;
1675	k = 0;
1676	for(i=0; i<=npoints-1; i++)
1677	{
1678	for(i_=0; i_<=nin-1;i_++)
1679	{
1680	network.x[i_] = xy[i,i_];
1681	}
1682	mlpprocess(ref network, ref network.x, ref network.y);
1683	if( mlpissoftmax(ref network) )
1684	{
1685
1686	//
1687	// class labels
1688	//
1689	lk = (int)Math.Round(xy[i,nin]);
1690	for(j=0; j<=nout-1; j++)
1691	{
1692	if( j==lk )
1693	{
1694	result = result+Math.Abs(1-network.y[j]);
1695	k = k+1;
1696	}
1697	}
1698	}
1699	else
1700	{
1701
1702	//
1703	// real outputs
1704	//
1705	for(j=0; j<=nout-1; j++)
1706	{
1707	if( (double)(xy[i,nin+j])!=(double)(0) )
1708	{
1709	result = result+Math.Abs(xy[i,nin+j]-network.y[j])/Math.Abs(xy[i,nin+j]);
1710	k = k+1;
1711	}
1712	}
1713	}
1714	}
1715	if( k!=0 )
1716	{
1717	result = result/k;
1718	}
1719	return result;
1720	}
1721
1722
1723	/*************************************************************************
1724	Gradient calculation. Internal subroutine.
1725
1726	-- ALGLIB --
1727	Copyright 04.11.2007 by Bochkanov Sergey
1728	*************************************************************************/
1729	public static void mlpgrad(ref multilayerperceptron network,
1730	ref double[] x,
1731	ref double[] desiredy,
1732	ref double e,
1733	ref double[] grad)
1734	{
1735	int i = 0;
1736	int nout = 0;
1737	int ntotal = 0;
1738
1739
1740	//
1741	// Prepare dError/dOut, internal structures
1742	//
1743	mlpprocess(ref network, ref x, ref network.y);
1744	nout = network.structinfo[2];
1745	ntotal = network.structinfo[3];
1746	e = 0;
1747	for(i=0; i<=ntotal-1; i++)
1748	{
1749	network.derror[i] = 0;
1750	}
1751	for(i=0; i<=nout-1; i++)
1752	{
1753	network.derror[ntotal-nout+i] = network.y[i]-desiredy[i];
1754	e = e+AP.Math.Sqr(network.y[i]-desiredy[i])/2;
1755	}
1756
1757	//
1758	// gradient
1759	//
1760	mlpinternalcalculategradient(ref network, ref network.neurons, ref network.weights, ref network.derror, ref grad, false);
1761	}
1762
1763
1764	/*************************************************************************
1765	Gradient calculation (natural error function). Internal subroutine.
1766
1767	-- ALGLIB --
1768	Copyright 04.11.2007 by Bochkanov Sergey
1769	*************************************************************************/
1770	public static void mlpgradn(ref multilayerperceptron network,
1771	ref double[] x,
1772	ref double[] desiredy,
1773	ref double e,
1774	ref double[] grad)
1775	{
1776	double s = 0;
1777	int i = 0;
1778	int nout = 0;
1779	int ntotal = 0;
1780
1781
1782	//
1783	// Prepare dError/dOut, internal structures
1784	//
1785	mlpprocess(ref network, ref x, ref network.y);
1786	nout = network.structinfo[2];
1787	ntotal = network.structinfo[3];
1788	for(i=0; i<=ntotal-1; i++)
1789	{
1790	network.derror[i] = 0;
1791	}
1792	e = 0;
1793	if( network.structinfo[6]==0 )
1794	{
1795
1796	//
1797	// Regression network, least squares
1798	//
1799	for(i=0; i<=nout-1; i++)
1800	{
1801	network.derror[ntotal-nout+i] = network.y[i]-desiredy[i];
1802	e = e+AP.Math.Sqr(network.y[i]-desiredy[i])/2;
1803	}
1804	}
1805	else
1806	{
1807
1808	//
1809	// Classification network, cross-entropy
1810	//
1811	s = 0;
1812	for(i=0; i<=nout-1; i++)
1813	{
1814	s = s+desiredy[i];
1815	}
1816	for(i=0; i<=nout-1; i++)
1817	{
1818	network.derror[ntotal-nout+i] = s*network.y[i]-desiredy[i];
1819	e = e+safecrossentropy(desiredy[i], network.y[i]);
1820	}
1821	}
1822
1823	//
1824	// gradient
1825	//
1826	mlpinternalcalculategradient(ref network, ref network.neurons, ref network.weights, ref network.derror, ref grad, true);
1827	}
1828
1829
1830	/*************************************************************************
1831	Batch gradient calculation. Internal subroutine.
1832
1833	-- ALGLIB --
1834	Copyright 04.11.2007 by Bochkanov Sergey
1835	*************************************************************************/
1836	public static void mlpgradbatch(ref multilayerperceptron network,
1837	ref double[,] xy,
1838	int ssize,
1839	ref double e,
1840	ref double[] grad)
1841	{
1842	int i = 0;
1843	int nin = 0;
1844	int nout = 0;
1845	int wcount = 0;
1846
1847	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1848	for(i=0; i<=wcount-1; i++)
1849	{
1850	grad[i] = 0;
1851	}
1852	e = 0;
1853	i = 0;
1854	while( i<=ssize-1 )
1855	{
1856	mlpchunkedgradient(ref network, ref xy, i, Math.Min(ssize, i+chunksize)-i, ref e, ref grad, false);
1857	i = i+chunksize;
1858	}
1859	}
1860
1861
1862	/*************************************************************************
1863	Batch gradient calculation (natural error function). Internal subroutine.
1864
1865	-- ALGLIB --
1866	Copyright 04.11.2007 by Bochkanov Sergey
1867	*************************************************************************/
1868	public static void mlpgradnbatch(ref multilayerperceptron network,
1869	ref double[,] xy,
1870	int ssize,
1871	ref double e,
1872	ref double[] grad)
1873	{
1874	int i = 0;
1875	int nin = 0;
1876	int nout = 0;
1877	int wcount = 0;
1878
1879	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1880	for(i=0; i<=wcount-1; i++)
1881	{
1882	grad[i] = 0;
1883	}
1884	e = 0;
1885	i = 0;
1886	while( i<=ssize-1 )
1887	{
1888	mlpchunkedgradient(ref network, ref xy, i, Math.Min(ssize, i+chunksize)-i, ref e, ref grad, true);
1889	i = i+chunksize;
1890	}
1891	}
1892
1893
1894	/*************************************************************************
1895	Batch Hessian calculation (natural error function) using R-algorithm.
1896	Internal subroutine.
1897
1898	-- ALGLIB --
1899	Copyright 26.01.2008 by Bochkanov Sergey.
1900
1901	Hessian calculation based on R-algorithm described in
1902	"Fast Exact Multiplication by the Hessian",
1903	B. A. Pearlmutter,
1904	Neural Computation, 1994.
1905	*************************************************************************/
1906	public static void mlphessiannbatch(ref multilayerperceptron network,
1907	ref double[,] xy,
1908	int ssize,
1909	ref double e,
1910	ref double[] grad,
1911	ref double[,] h)
1912	{
1913	mlphessianbatchinternal(ref network, ref xy, ssize, true, ref e, ref grad, ref h);
1914	}
1915
1916
1917	/*************************************************************************
1918	Batch Hessian calculation using R-algorithm.
1919	Internal subroutine.
1920
1921	-- ALGLIB --
1922	Copyright 26.01.2008 by Bochkanov Sergey.
1923
1924	Hessian calculation based on R-algorithm described in
1925	"Fast Exact Multiplication by the Hessian",
1926	B. A. Pearlmutter,
1927	Neural Computation, 1994.
1928	*************************************************************************/
1929	public static void mlphessianbatch(ref multilayerperceptron network,
1930	ref double[,] xy,
1931	int ssize,
1932	ref double e,
1933	ref double[] grad,
1934	ref double[,] h)
1935	{
1936	mlphessianbatchinternal(ref network, ref xy, ssize, false, ref e, ref grad, ref h);
1937	}
1938
1939
1940	/*************************************************************************
1941	Internal subroutine, shouldn't be called by user.
1942	*************************************************************************/
1943	public static void mlpinternalprocessvector(ref int[] structinfo,
1944	ref double[] weights,
1945	ref double[] columnmeans,
1946	ref double[] columnsigmas,
1947	ref double[] neurons,
1948	ref double[] dfdnet,
1949	ref double[] x,
1950	ref double[] y)
1951	{
1952	int i = 0;
1953	int j = 0;
1954	int n1 = 0;
1955	int n2 = 0;
1956	int w1 = 0;
1957	int w2 = 0;
1958	int ntotal = 0;
1959	int nin = 0;
1960	int nout = 0;
1961	int istart = 0;
1962	int offs = 0;
1963	double net = 0;
1964	double e = 0;
1965	double f = 0;
1966	double df = 0;
1967	double d2f = 0;
1968	double mx = 0;
1969	bool perr = new bool();
1970	int i_ = 0;
1971	int i1_ = 0;
1972
1973
1974	//
1975	// Read network geometry
1976	//
1977	nin = structinfo[1];
1978	nout = structinfo[2];
1979	ntotal = structinfo[3];
1980	istart = structinfo[5];
1981
1982	//
1983	// Inputs standartisation and putting in the network
1984	//
1985	for(i=0; i<=nin-1; i++)
1986	{
1987	if( (double)(columnsigmas[i])!=(double)(0) )
1988	{
1989	neurons[i] = (x[i]-columnmeans[i])/columnsigmas[i];
1990	}
1991	else
1992	{
1993	neurons[i] = x[i]-columnmeans[i];
1994	}
1995	}
1996
1997	//
1998	// Process network
1999	//
2000	for(i=0; i<=ntotal-1; i++)
2001	{
2002	offs = istart+i*nfieldwidth;
2003	if( structinfo[offs+0]>0 )
2004	{
2005
2006	//
2007	// Activation function
2008	//
2009	mlpactivationfunction(neurons[structinfo[offs+2]], structinfo[offs+0], ref f, ref df, ref d2f);
2010	neurons[i] = f;
2011	dfdnet[i] = df;
2012	}
2013	if( structinfo[offs+0]==0 )
2014	{
2015
2016	//
2017	// Adaptive summator
2018	//
2019	n1 = structinfo[offs+2];
2020	n2 = n1+structinfo[offs+1]-1;
2021	w1 = structinfo[offs+3];
2022	w2 = w1+structinfo[offs+1]-1;
2023	i1_ = (n1)-(w1);
2024	net = 0.0;
2025	for(i_=w1; i_<=w2;i_++)
2026	{
2027	net += weights[i_]*neurons[i_+i1_];
2028	}
2029	neurons[i] = net;
2030	dfdnet[i] = 1.0;
2031	}
2032	if( structinfo[offs+0]<0 )
2033	{
2034	perr = true;
2035	if( structinfo[offs+0]==-2 )
2036	{
2037
2038	//
2039	// input neuron, left unchanged
2040	//
2041	perr = false;
2042	}
2043	if( structinfo[offs+0]==-3 )
2044	{
2045
2046	//
2047	// "-1" neuron
2048	//
2049	neurons[i] = -1;
2050	perr = false;
2051	}
2052	if( structinfo[offs+0]==-4 )
2053	{
2054
2055	//
2056	// "0" neuron
2057	//
2058	neurons[i] = 0;
2059	perr = false;
2060	}
2061	System.Diagnostics.Debug.Assert(!perr, "MLPInternalProcessVector: internal error - unknown neuron type!");
2062	}
2063	}
2064
2065	//
2066	// Extract result
2067	//
2068	i1_ = (ntotal-nout) - (0);
2069	for(i_=0; i_<=nout-1;i_++)
2070	{
2071	y[i_] = neurons[i_+i1_];
2072	}
2073
2074	//
2075	// Softmax post-processing or standardisation if needed
2076	//
2077	System.Diagnostics.Debug.Assert(structinfo[6]==0 \| structinfo[6]==1, "MLPInternalProcessVector: unknown normalization type!");
2078	if( structinfo[6]==1 )
2079	{
2080
2081	//
2082	// Softmax
2083	//
2084	mx = y[0];
2085	for(i=1; i<=nout-1; i++)
2086	{
2087	mx = Math.Max(mx, y[i]);
2088	}
2089	net = 0;
2090	for(i=0; i<=nout-1; i++)
2091	{
2092	y[i] = Math.Exp(y[i]-mx);
2093	net = net+y[i];
2094	}
2095	for(i=0; i<=nout-1; i++)
2096	{
2097	y[i] = y[i]/net;
2098	}
2099	}
2100	else
2101	{
2102
2103	//
2104	// Standardisation
2105	//
2106	for(i=0; i<=nout-1; i++)
2107	{
2108	y[i] = y[i]*columnsigmas[nin+i]+columnmeans[nin+i];
2109	}
2110	}
2111	}
2112
2113
2114	/*************************************************************************
2115	Internal subroutine: adding new input layer to network
2116	*************************************************************************/
2117	private static void addinputlayer(int ncount,
2118	ref int[] lsizes,
2119	ref int[] ltypes,
2120	ref int[] lconnfirst,
2121	ref int[] lconnlast,
2122	ref int lastproc)
2123	{
2124	lsizes[0] = ncount;
2125	ltypes[0] = -2;
2126	lconnfirst[0] = 0;
2127	lconnlast[0] = 0;
2128	lastproc = 0;
2129	}
2130
2131
2132	/*************************************************************************
2133	Internal subroutine: adding new summator layer to network
2134	*************************************************************************/
2135	private static void addbiasedsummatorlayer(int ncount,
2136	ref int[] lsizes,
2137	ref int[] ltypes,
2138	ref int[] lconnfirst,
2139	ref int[] lconnlast,
2140	ref int lastproc)
2141	{
2142	lsizes[lastproc+1] = 1;
2143	ltypes[lastproc+1] = -3;
2144	lconnfirst[lastproc+1] = 0;
2145	lconnlast[lastproc+1] = 0;
2146	lsizes[lastproc+2] = ncount;
2147	ltypes[lastproc+2] = 0;
2148	lconnfirst[lastproc+2] = lastproc;
2149	lconnlast[lastproc+2] = lastproc+1;
2150	lastproc = lastproc+2;
2151	}
2152
2153
2154	/*************************************************************************
2155	Internal subroutine: adding new summator layer to network
2156	*************************************************************************/
2157	private static void addactivationlayer(int functype,
2158	ref int[] lsizes,
2159	ref int[] ltypes,
2160	ref int[] lconnfirst,
2161	ref int[] lconnlast,
2162	ref int lastproc)
2163	{
2164	System.Diagnostics.Debug.Assert(functype>0, "AddActivationLayer: incorrect function type");
2165	lsizes[lastproc+1] = lsizes[lastproc];
2166	ltypes[lastproc+1] = functype;
2167	lconnfirst[lastproc+1] = lastproc;
2168	lconnlast[lastproc+1] = lastproc;
2169	lastproc = lastproc+1;
2170	}
2171
2172
2173	/*************************************************************************
2174	Internal subroutine: adding new zero layer to network
2175	*************************************************************************/
2176	private static void addzerolayer(ref int[] lsizes,
2177	ref int[] ltypes,
2178	ref int[] lconnfirst,
2179	ref int[] lconnlast,
2180	ref int lastproc)
2181	{
2182	lsizes[lastproc+1] = 1;
2183	ltypes[lastproc+1] = -4;
2184	lconnfirst[lastproc+1] = 0;
2185	lconnlast[lastproc+1] = 0;
2186	lastproc = lastproc+1;
2187	}
2188
2189
2190	/*************************************************************************
2191	Internal subroutine.
2192
2193	-- ALGLIB --
2194	Copyright 04.11.2007 by Bochkanov Sergey
2195	*************************************************************************/
2196	private static void mlpcreate(int nin,
2197	int nout,
2198	ref int[] lsizes,
2199	ref int[] ltypes,
2200	ref int[] lconnfirst,
2201	ref int[] lconnlast,
2202	int layerscount,
2203	bool isclsnet,
2204	ref multilayerperceptron network)
2205	{
2206	int i = 0;
2207	int j = 0;
2208	int ssize = 0;
2209	int ntotal = 0;
2210	int wcount = 0;
2211	int offs = 0;
2212	int nprocessed = 0;
2213	int wallocated = 0;
2214	int[] localtemp = new int[0];
2215	int[] lnfirst = new int[0];
2216	int[] lnsyn = new int[0];
2217
2218
2219	//
2220	// Check
2221	//
2222	System.Diagnostics.Debug.Assert(layerscount>0, "MLPCreate: wrong parameters!");
2223	System.Diagnostics.Debug.Assert(ltypes[0]==-2, "MLPCreate: wrong LTypes[0] (must be -2)!");
2224	for(i=0; i<=layerscount-1; i++)
2225	{
2226	System.Diagnostics.Debug.Assert(lsizes[i]>0, "MLPCreate: wrong LSizes!");
2227	System.Diagnostics.Debug.Assert(lconnfirst[i]>=0 & (lconnfirst[i]<i \| i==0), "MLPCreate: wrong LConnFirst!");
2228	System.Diagnostics.Debug.Assert(lconnlast[i]>=lconnfirst[i] & (lconnlast[i]<i \| i==0), "MLPCreate: wrong LConnLast!");
2229	}
2230
2231	//
2232	// Build network geometry
2233	//
2234	lnfirst = new int[layerscount-1+1];
2235	lnsyn = new int[layerscount-1+1];
2236	ntotal = 0;
2237	wcount = 0;
2238	for(i=0; i<=layerscount-1; i++)
2239	{
2240
2241	//
2242	// Analyze connections.
2243	// This code must throw an assertion in case of unknown LTypes[I]
2244	//
2245	lnsyn[i] = -1;
2246	if( ltypes[i]>=0 )
2247	{
2248	lnsyn[i] = 0;
2249	for(j=lconnfirst[i]; j<=lconnlast[i]; j++)
2250	{
2251	lnsyn[i] = lnsyn[i]+lsizes[j];
2252	}
2253	}
2254	else
2255	{
2256	if( ltypes[i]==-2 \| ltypes[i]==-3 \| ltypes[i]==-4 )
2257	{
2258	lnsyn[i] = 0;
2259	}
2260	}
2261	System.Diagnostics.Debug.Assert(lnsyn[i]>=0, "MLPCreate: internal error #0!");
2262
2263	//
2264	// Other info
2265	//
2266	lnfirst[i] = ntotal;
2267	ntotal = ntotal+lsizes[i];
2268	if( ltypes[i]==0 )
2269	{
2270	wcount = wcount+lnsyn[i]*lsizes[i];
2271	}
2272	}
2273	ssize = 7+ntotal*nfieldwidth;
2274
2275	//
2276	// Allocate
2277	//
2278	network.structinfo = new int[ssize-1+1];
2279	network.weights = new double[wcount-1+1];
2280	if( isclsnet )
2281	{
2282	network.columnmeans = new double[nin-1+1];
2283	network.columnsigmas = new double[nin-1+1];
2284	}
2285	else
2286	{
2287	network.columnmeans = new double[nin+nout-1+1];
2288	network.columnsigmas = new double[nin+nout-1+1];
2289	}
2290	network.neurons = new double[ntotal-1+1];
2291	network.chunks = new double[3*ntotal+1, chunksize-1+1];
2292	network.nwbuf = new double[Math.Max(wcount, 2*nout)-1+1];
2293	network.dfdnet = new double[ntotal-1+1];
2294	network.x = new double[nin-1+1];
2295	network.y = new double[nout-1+1];
2296	network.derror = new double[ntotal-1+1];
2297
2298	//
2299	// Fill structure: global info
2300	//
2301	network.structinfo[0] = ssize;
2302	network.structinfo[1] = nin;
2303	network.structinfo[2] = nout;
2304	network.structinfo[3] = ntotal;
2305	network.structinfo[4] = wcount;
2306	network.structinfo[5] = 7;
2307	if( isclsnet )
2308	{
2309	network.structinfo[6] = 1;
2310	}
2311	else
2312	{
2313	network.structinfo[6] = 0;
2314	}
2315
2316	//
2317	// Fill structure: neuron connections
2318	//
2319	nprocessed = 0;
2320	wallocated = 0;
2321	for(i=0; i<=layerscount-1; i++)
2322	{
2323	for(j=0; j<=lsizes[i]-1; j++)
2324	{
2325	offs = network.structinfo[5]+nprocessed*nfieldwidth;
2326	network.structinfo[offs+0] = ltypes[i];
2327	if( ltypes[i]==0 )
2328	{
2329
2330	//
2331	// Adaptive summator:
2332	// * connections with weights to previous neurons
2333	//
2334	network.structinfo[offs+1] = lnsyn[i];
2335	network.structinfo[offs+2] = lnfirst[lconnfirst[i]];
2336	network.structinfo[offs+3] = wallocated;
2337	wallocated = wallocated+lnsyn[i];
2338	nprocessed = nprocessed+1;
2339	}
2340	if( ltypes[i]>0 )
2341	{
2342
2343	//
2344	// Activation layer:
2345	// * each neuron connected to one (only one) of previous neurons.
2346	// * no weights
2347	//
2348	network.structinfo[offs+1] = 1;
2349	network.structinfo[offs+2] = lnfirst[lconnfirst[i]]+j;
2350	network.structinfo[offs+3] = -1;
2351	nprocessed = nprocessed+1;
2352	}
2353	if( ltypes[i]==-2 \| ltypes[i]==-3 \| ltypes[i]==-4 )
2354	{
2355	nprocessed = nprocessed+1;
2356	}
2357	}
2358	}
2359	System.Diagnostics.Debug.Assert(wallocated==wcount, "MLPCreate: internal error #1!");
2360	System.Diagnostics.Debug.Assert(nprocessed==ntotal, "MLPCreate: internal error #2!");
2361
2362	//
2363	// Fill weights by small random values
2364	// Initialize means and sigmas
2365	//
2366	for(i=0; i<=wcount-1; i++)
2367	{
2368	network.weights[i] = AP.Math.RandomReal()-0.5;
2369	}
2370	for(i=0; i<=nin-1; i++)
2371	{
2372	network.columnmeans[i] = 0;
2373	network.columnsigmas[i] = 1;
2374	}
2375	if( !isclsnet )
2376	{
2377	for(i=0; i<=nout-1; i++)
2378	{
2379	network.columnmeans[nin+i] = 0;
2380	network.columnsigmas[nin+i] = 1;
2381	}
2382	}
2383	}
2384
2385
2386	/*************************************************************************
2387	Internal subroutine
2388
2389	-- ALGLIB --
2390	Copyright 04.11.2007 by Bochkanov Sergey
2391	*************************************************************************/
2392	private static void mlpactivationfunction(double net,
2393	int k,
2394	ref double f,
2395	ref double df,
2396	ref double d2f)
2397	{
2398	double net2 = 0;
2399	double arg = 0;
2400	double root = 0;
2401	double r = 0;
2402
2403	f = 0;
2404	df = 0;
2405	if( k==1 )
2406	{
2407
2408	//
2409	// TanH activation function
2410	//
2411	if( (double)(Math.Abs(net))<(double)(100) )
2412	{
2413	f = Math.Tanh(net);
2414	}
2415	else
2416	{
2417	f = Math.Sign(net);
2418	}
2419	df = 1-AP.Math.Sqr(f);
2420	d2f = -(2fdf);
2421	return;
2422	}
2423	if( k==3 )
2424	{
2425
2426	//
2427	// EX activation function
2428	//
2429	if( (double)(net)>=(double)(0) )
2430	{
2431	net2 = net*net;
2432	arg = net2+1;
2433	root = Math.Sqrt(arg);
2434	f = net+root;
2435	r = net/root;
2436	df = 1+r;
2437	d2f = (root-net*r)/arg;
2438	}
2439	else
2440	{
2441	f = Math.Exp(net);
2442	df = f;
2443	d2f = f;
2444	}
2445	return;
2446	}
2447	if( k==2 )
2448	{
2449	f = Math.Exp(-AP.Math.Sqr(net));
2450	df = -(2netf);
2451	d2f = -(2(f+dfnet));
2452	return;
2453	}
2454	}
2455
2456
2457	/*************************************************************************
2458	Internal subroutine for Hessian calculation.
2459
2460	WARNING!!! Unspeakable math far beyong human capabilities :)
2461	*************************************************************************/
2462	private static void mlphessianbatchinternal(ref multilayerperceptron network,
2463	ref double[,] xy,
2464	int ssize,
2465	bool naturalerr,
2466	ref double e,
2467	ref double[] grad,
2468	ref double[,] h)
2469	{
2470	int nin = 0;
2471	int nout = 0;
2472	int wcount = 0;
2473	int ntotal = 0;
2474	int istart = 0;
2475	int i = 0;
2476	int j = 0;
2477	int k = 0;
2478	int kl = 0;
2479	int offs = 0;
2480	int n1 = 0;
2481	int n2 = 0;
2482	int w1 = 0;
2483	int w2 = 0;
2484	double s = 0;
2485	double t = 0;
2486	double v = 0;
2487	double et = 0;
2488	bool bflag = new bool();
2489	double f = 0;
2490	double df = 0;
2491	double d2f = 0;
2492	double deidyj = 0;
2493	double mx = 0;
2494	double net = 0;
2495	double q = 0;
2496	double z = 0;
2497	double s2 = 0;
2498	double expi = 0;
2499	double expj = 0;
2500	double[] x = new double[0];
2501	double[] desiredy = new double[0];
2502	double[] gt = new double[0];
2503	double[] zeros = new double[0];
2504	double[,] rx = new double[0,0];
2505	double[,] ry = new double[0,0];
2506	double[,] rdx = new double[0,0];
2507	double[,] rdy = new double[0,0];
2508	int i_ = 0;
2509	int i1_ = 0;
2510
2511	mlpproperties(ref network, ref nin, ref nout, ref wcount);
2512	ntotal = network.structinfo[3];
2513	istart = network.structinfo[5];
2514
2515	//
2516	// Prepare
2517	//
2518	x = new double[nin-1+1];
2519	desiredy = new double[nout-1+1];
2520	zeros = new double[wcount-1+1];
2521	gt = new double[wcount-1+1];
2522	rx = new double[ntotal+nout-1+1, wcount-1+1];
2523	ry = new double[ntotal+nout-1+1, wcount-1+1];
2524	rdx = new double[ntotal+nout-1+1, wcount-1+1];
2525	rdy = new double[ntotal+nout-1+1, wcount-1+1];
2526	e = 0;
2527	for(i=0; i<=wcount-1; i++)
2528	{
2529	zeros[i] = 0;
2530	}
2531	for(i_=0; i_<=wcount-1;i_++)
2532	{
2533	grad[i_] = zeros[i_];
2534	}
2535	for(i=0; i<=wcount-1; i++)
2536	{
2537	for(i_=0; i_<=wcount-1;i_++)
2538	{
2539	h[i,i_] = zeros[i_];
2540	}
2541	}
2542
2543	//
2544	// Process
2545	//
2546	for(k=0; k<=ssize-1; k++)
2547	{
2548
2549	//
2550	// Process vector with MLPGradN.
2551	// Now Neurons, DFDNET and DError contains results of the last run.
2552	//
2553	for(i_=0; i_<=nin-1;i_++)
2554	{
2555	x[i_] = xy[k,i_];
2556	}
2557	if( mlpissoftmax(ref network) )
2558	{
2559
2560	//
2561	// class labels outputs
2562	//
2563	kl = (int)Math.Round(xy[k,nin]);
2564	for(i=0; i<=nout-1; i++)
2565	{
2566	if( i==kl )
2567	{
2568	desiredy[i] = 1;
2569	}
2570	else
2571	{
2572	desiredy[i] = 0;
2573	}
2574	}
2575	}
2576	else
2577	{
2578
2579	//
2580	// real outputs
2581	//
2582	i1_ = (nin) - (0);
2583	for(i_=0; i_<=nout-1;i_++)
2584	{
2585	desiredy[i_] = xy[k,i_+i1_];
2586	}
2587	}
2588	if( naturalerr )
2589	{
2590	mlpgradn(ref network, ref x, ref desiredy, ref et, ref gt);
2591	}
2592	else
2593	{
2594	mlpgrad(ref network, ref x, ref desiredy, ref et, ref gt);
2595	}
2596
2597	//
2598	// grad, error
2599	//
2600	e = e+et;
2601	for(i_=0; i_<=wcount-1;i_++)
2602	{
2603	grad[i_] = grad[i_] + gt[i_];
2604	}
2605
2606	//
2607	// Hessian.
2608	// Forward pass of the R-algorithm
2609	//
2610	for(i=0; i<=ntotal-1; i++)
2611	{
2612	offs = istart+i*nfieldwidth;
2613	for(i_=0; i_<=wcount-1;i_++)
2614	{
2615	rx[i,i_] = zeros[i_];
2616	}
2617	for(i_=0; i_<=wcount-1;i_++)
2618	{
2619	ry[i,i_] = zeros[i_];
2620	}
2621	if( network.structinfo[offs+0]>0 )
2622	{
2623
2624	//
2625	// Activation function
2626	//
2627	n1 = network.structinfo[offs+2];
2628	for(i_=0; i_<=wcount-1;i_++)
2629	{
2630	rx[i,i_] = ry[n1,i_];
2631	}
2632	v = network.dfdnet[i];
2633	for(i_=0; i_<=wcount-1;i_++)
2634	{
2635	ry[i,i_] = v*rx[i,i_];
2636	}
2637	}
2638	if( network.structinfo[offs+0]==0 )
2639	{
2640
2641	//
2642	// Adaptive summator
2643	//
2644	n1 = network.structinfo[offs+2];
2645	n2 = n1+network.structinfo[offs+1]-1;
2646	w1 = network.structinfo[offs+3];
2647	w2 = w1+network.structinfo[offs+1]-1;
2648	for(j=n1; j<=n2; j++)
2649	{
2650	v = network.weights[w1+j-n1];
2651	for(i_=0; i_<=wcount-1;i_++)
2652	{
2653	rx[i,i_] = rx[i,i_] + v*ry[j,i_];
2654	}
2655	rx[i,w1+j-n1] = rx[i,w1+j-n1]+network.neurons[j];
2656	}
2657	for(i_=0; i_<=wcount-1;i_++)
2658	{
2659	ry[i,i_] = rx[i,i_];
2660	}
2661	}
2662	if( network.structinfo[offs+0]<0 )
2663	{
2664	bflag = true;
2665	if( network.structinfo[offs+0]==-2 )
2666	{
2667
2668	//
2669	// input neuron, left unchanged
2670	//
2671	bflag = false;
2672	}
2673	if( network.structinfo[offs+0]==-3 )
2674	{
2675
2676	//
2677	// "-1" neuron, left unchanged
2678	//
2679	bflag = false;
2680	}
2681	if( network.structinfo[offs+0]==-4 )
2682	{
2683
2684	//
2685	// "0" neuron, left unchanged
2686	//
2687	bflag = false;
2688	}
2689	System.Diagnostics.Debug.Assert(!bflag, "MLPHessianNBatch: internal error - unknown neuron type!");
2690	}
2691	}
2692
2693	//
2694	// Hessian. Backward pass of the R-algorithm.
2695	//
2696	// Stage 1. Initialize RDY
2697	//
2698	for(i=0; i<=ntotal+nout-1; i++)
2699	{
2700	for(i_=0; i_<=wcount-1;i_++)
2701	{
2702	rdy[i,i_] = zeros[i_];
2703	}
2704	}
2705	if( network.structinfo[6]==0 )
2706	{
2707
2708	//
2709	// Standardisation.
2710	//
2711	// In context of the Hessian calculation standardisation
2712	// is considered as additional layer with weightless
2713	// activation function:
2714	//
2715	// F(NET) := Sigma*NET
2716	//
2717	// So we add one more layer to forward pass, and
2718	// make forward/backward pass through this layer.
2719	//
2720	for(i=0; i<=nout-1; i++)
2721	{
2722	n1 = ntotal-nout+i;
2723	n2 = ntotal+i;
2724
2725	//
2726	// Forward pass from N1 to N2
2727	//
2728	for(i_=0; i_<=wcount-1;i_++)
2729	{
2730	rx[n2,i_] = ry[n1,i_];
2731	}
2732	v = network.columnsigmas[nin+i];
2733	for(i_=0; i_<=wcount-1;i_++)
2734	{
2735	ry[n2,i_] = v*rx[n2,i_];
2736	}
2737
2738	//
2739	// Initialization of RDY
2740	//
2741	for(i_=0; i_<=wcount-1;i_++)
2742	{
2743	rdy[n2,i_] = ry[n2,i_];
2744	}
2745
2746	//
2747	// Backward pass from N2 to N1:
2748	// 1. Calculate R(dE/dX).
2749	// 2. No R(dE/dWij) is needed since weight of activation neuron
2750	// is fixed to 1. So we can update R(dE/dY) for
2751	// the connected neuron (note that Vij=0, Wij=1)
2752	//
2753	df = network.columnsigmas[nin+i];
2754	for(i_=0; i_<=wcount-1;i_++)
2755	{
2756	rdx[n2,i_] = df*rdy[n2,i_];
2757	}
2758	for(i_=0; i_<=wcount-1;i_++)
2759	{
2760	rdy[n1,i_] = rdy[n1,i_] + rdx[n2,i_];
2761	}
2762	}
2763	}
2764	else
2765	{
2766
2767	//
2768	// Softmax.
2769	//
2770	// Initialize RDY using generalized expression for ei'(yi)
2771	// (see expression (9) from p. 5 of "Fast Exact Multiplication by the Hessian").
2772	//
2773	// When we are working with softmax network, generalized
2774	// expression for ei'(yi) is used because softmax
2775	// normalization leads to ei, which depends on all y's
2776	//
2777	if( naturalerr )
2778	{
2779
2780	//
2781	// softmax + cross-entropy.
2782	// We have:
2783	//
2784	// S = sum(exp(yk)),
2785	// ei = sum(trn)*exp(yi)/S-trn_i
2786	//
2787	// j=i: d(ei)/d(yj) = Texp(yi)(S-exp(yi))/S^2
2788	// j<>i: d(ei)/d(yj) = -Texp(yi)exp(yj)/S^2
2789	//
2790	t = 0;
2791	for(i=0; i<=nout-1; i++)
2792	{
2793	t = t+desiredy[i];
2794	}
2795	mx = network.neurons[ntotal-nout];
2796	for(i=0; i<=nout-1; i++)
2797	{
2798	mx = Math.Max(mx, network.neurons[ntotal-nout+i]);
2799	}
2800	s = 0;
2801	for(i=0; i<=nout-1; i++)
2802	{
2803	network.nwbuf[i] = Math.Exp(network.neurons[ntotal-nout+i]-mx);
2804	s = s+network.nwbuf[i];
2805	}
2806	for(i=0; i<=nout-1; i++)
2807	{
2808	for(j=0; j<=nout-1; j++)
2809	{
2810	if( j==i )
2811	{
2812	deidyj = tnetwork.nwbuf[i](s-network.nwbuf[i])/AP.Math.Sqr(s);
2813	for(i_=0; i_<=wcount-1;i_++)
2814	{
2815	rdy[ntotal-nout+i,i_] = rdy[ntotal-nout+i,i_] + deidyj*ry[ntotal-nout+i,i_];
2816	}
2817	}
2818	else
2819	{
2820	deidyj = -(tnetwork.nwbuf[i]network.nwbuf[j]/AP.Math.Sqr(s));
2821	for(i_=0; i_<=wcount-1;i_++)
2822	{
2823	rdy[ntotal-nout+i,i_] = rdy[ntotal-nout+i,i_] + deidyj*ry[ntotal-nout+j,i_];
2824	}
2825	}
2826	}
2827	}
2828	}
2829	else
2830	{
2831
2832	//
2833	// For a softmax + squared error we have expression
2834	// far beyond human imagination so we dont even try
2835	// to comment on it. Just enjoy the code...
2836	//
2837	// P.S. That's why "natural error" is called "natural" -
2838	// compact beatiful expressions, fast code....
2839	//
2840	mx = network.neurons[ntotal-nout];
2841	for(i=0; i<=nout-1; i++)
2842	{
2843	mx = Math.Max(mx, network.neurons[ntotal-nout+i]);
2844	}
2845	s = 0;
2846	s2 = 0;
2847	for(i=0; i<=nout-1; i++)
2848	{
2849	network.nwbuf[i] = Math.Exp(network.neurons[ntotal-nout+i]-mx);
2850	s = s+network.nwbuf[i];
2851	s2 = s2+AP.Math.Sqr(network.nwbuf[i]);
2852	}
2853	q = 0;
2854	for(i=0; i<=nout-1; i++)
2855	{
2856	q = q+(network.y[i]-desiredy[i])*network.nwbuf[i];
2857	}
2858	for(i=0; i<=nout-1; i++)
2859	{
2860	z = -q+(network.y[i]-desiredy[i])*s;
2861	expi = network.nwbuf[i];
2862	for(j=0; j<=nout-1; j++)
2863	{
2864	expj = network.nwbuf[j];
2865	if( j==i )
2866	{
2867	deidyj = expi/AP.Math.Sqr(s)((z+expi)(s-2expi)/s+expis2/AP.Math.Sqr(s));
2868	}
2869	else
2870	{
2871	deidyj = expiexpj/AP.Math.Sqr(s)(s2/AP.Math.Sqr(s)-2*z/s-(expi+expj)/s+(network.y[i]-desiredy[i])-(network.y[j]-desiredy[j]));
2872	}
2873	for(i_=0; i_<=wcount-1;i_++)
2874	{
2875	rdy[ntotal-nout+i,i_] = rdy[ntotal-nout+i,i_] + deidyj*ry[ntotal-nout+j,i_];
2876	}
2877	}
2878	}
2879	}
2880	}
2881
2882	//
2883	// Hessian. Backward pass of the R-algorithm
2884	//
2885	// Stage 2. Process.
2886	//
2887	for(i=ntotal-1; i>=0; i--)
2888	{
2889
2890	//
2891	// Possible variants:
2892	// 1. Activation function
2893	// 2. Adaptive summator
2894	// 3. Special neuron
2895	//
2896	offs = istart+i*nfieldwidth;
2897	if( network.structinfo[offs+0]>0 )
2898	{
2899	n1 = network.structinfo[offs+2];
2900
2901	//
2902	// First, calculate R(dE/dX).
2903	//
2904	mlpactivationfunction(network.neurons[n1], network.structinfo[offs+0], ref f, ref df, ref d2f);
2905	v = d2f*network.derror[i];
2906	for(i_=0; i_<=wcount-1;i_++)
2907	{
2908	rdx[i,i_] = df*rdy[i,i_];
2909	}
2910	for(i_=0; i_<=wcount-1;i_++)
2911	{
2912	rdx[i,i_] = rdx[i,i_] + v*rx[i,i_];
2913	}
2914
2915	//
2916	// No R(dE/dWij) is needed since weight of activation neuron
2917	// is fixed to 1.
2918	//
2919	// So we can update R(dE/dY) for the connected neuron.
2920	// (note that Vij=0, Wij=1)
2921	//
2922	for(i_=0; i_<=wcount-1;i_++)
2923	{
2924	rdy[n1,i_] = rdy[n1,i_] + rdx[i,i_];
2925	}
2926	}
2927	if( network.structinfo[offs+0]==0 )
2928	{
2929
2930	//
2931	// Adaptive summator
2932	//
2933	n1 = network.structinfo[offs+2];
2934	n2 = n1+network.structinfo[offs+1]-1;
2935	w1 = network.structinfo[offs+3];
2936	w2 = w1+network.structinfo[offs+1]-1;
2937
2938	//
2939	// First, calculate R(dE/dX).
2940	//
2941	for(i_=0; i_<=wcount-1;i_++)
2942	{
2943	rdx[i,i_] = rdy[i,i_];
2944	}
2945
2946	//
2947	// Then, calculate R(dE/dWij)
2948	//
2949	for(j=w1; j<=w2; j++)
2950	{
2951	v = network.neurons[n1+j-w1];
2952	for(i_=0; i_<=wcount-1;i_++)
2953	{
2954	h[j,i_] = h[j,i_] + v*rdx[i,i_];
2955	}
2956	v = network.derror[i];
2957	for(i_=0; i_<=wcount-1;i_++)
2958	{
2959	h[j,i_] = h[j,i_] + v*ry[n1+j-w1,i_];
2960	}
2961	}
2962
2963	//
2964	// And finally, update R(dE/dY) for connected neurons.
2965	//
2966	for(j=w1; j<=w2; j++)
2967	{
2968	v = network.weights[j];
2969	for(i_=0; i_<=wcount-1;i_++)
2970	{
2971	rdy[n1+j-w1,i_] = rdy[n1+j-w1,i_] + v*rdx[i,i_];
2972	}
2973	rdy[n1+j-w1,j] = rdy[n1+j-w1,j]+network.derror[i];
2974	}
2975	}
2976	if( network.structinfo[offs+0]<0 )
2977	{
2978	bflag = false;
2979	if( network.structinfo[offs+0]==-2 \| network.structinfo[offs+0]==-3 \| network.structinfo[offs+0]==-4 )
2980	{
2981
2982	//
2983	// Special neuron type, no back-propagation required
2984	//
2985	bflag = true;
2986	}
2987	System.Diagnostics.Debug.Assert(bflag, "MLPHessianNBatch: unknown neuron type!");
2988	}
2989	}
2990	}
2991	}
2992
2993
2994	/*************************************************************************
2995	Internal subroutine
2996
2997	Network must be processed by MLPProcess on X
2998	*************************************************************************/
2999	private static void mlpinternalcalculategradient(ref multilayerperceptron network,
3000	ref double[] neurons,
3001	ref double[] weights,
3002	ref double[] derror,
3003	ref double[] grad,
3004	bool naturalerrorfunc)
3005	{
3006	int i = 0;
3007	int j = 0;
3008	int n1 = 0;
3009	int n2 = 0;
3010	int w1 = 0;
3011	int w2 = 0;
3012	int ntotal = 0;
3013	int istart = 0;
3014	int nin = 0;
3015	int nout = 0;
3016	int offs = 0;
3017	double dedf = 0;
3018	double dfdnet = 0;
3019	double v = 0;
3020	double fown = 0;
3021	double deown = 0;
3022	double net = 0;
3023	double mx = 0;
3024	bool bflag = new bool();
3025	int i_ = 0;
3026	int i1_ = 0;
3027
3028
3029	//
3030	// Read network geometry
3031	//
3032	nin = network.structinfo[1];
3033	nout = network.structinfo[2];
3034	ntotal = network.structinfo[3];
3035	istart = network.structinfo[5];
3036
3037	//
3038	// Pre-processing of dError/dOut:
3039	// from dError/dOut(normalized) to dError/dOut(non-normalized)
3040	//
3041	System.Diagnostics.Debug.Assert(network.structinfo[6]==0 \| network.structinfo[6]==1, "MLPInternalCalculateGradient: unknown normalization type!");
3042	if( network.structinfo[6]==1 )
3043	{
3044
3045	//
3046	// Softmax
3047	//
3048	if( !naturalerrorfunc )
3049	{
3050	mx = network.neurons[ntotal-nout];
3051	for(i=0; i<=nout-1; i++)
3052	{
3053	mx = Math.Max(mx, network.neurons[ntotal-nout+i]);
3054	}
3055	net = 0;
3056	for(i=0; i<=nout-1; i++)
3057	{
3058	network.nwbuf[i] = Math.Exp(network.neurons[ntotal-nout+i]-mx);
3059	net = net+network.nwbuf[i];
3060	}
3061	i1_ = (0)-(ntotal-nout);
3062	v = 0.0;
3063	for(i_=ntotal-nout; i_<=ntotal-1;i_++)
3064	{
3065	v += network.derror[i_]*network.nwbuf[i_+i1_];
3066	}
3067	for(i=0; i<=nout-1; i++)
3068	{
3069	fown = network.nwbuf[i];
3070	deown = network.derror[ntotal-nout+i];
3071	network.nwbuf[nout+i] = (-v+deownfown+deown(net-fown))*fown/AP.Math.Sqr(net);
3072	}
3073	for(i=0; i<=nout-1; i++)
3074	{
3075	network.derror[ntotal-nout+i] = network.nwbuf[nout+i];
3076	}
3077	}
3078	}
3079	else
3080	{
3081
3082	//
3083	// Un-standardisation
3084	//
3085	for(i=0; i<=nout-1; i++)
3086	{
3087	network.derror[ntotal-nout+i] = network.derror[ntotal-nout+i]*network.columnsigmas[nin+i];
3088	}
3089	}
3090
3091	//
3092	// Backpropagation
3093	//
3094	for(i=ntotal-1; i>=0; i--)
3095	{
3096
3097	//
3098	// Extract info
3099	//
3100	offs = istart+i*nfieldwidth;
3101	if( network.structinfo[offs+0]>0 )
3102	{
3103
3104	//
3105	// Activation function
3106	//
3107	dedf = network.derror[i];
3108	dfdnet = network.dfdnet[i];
3109	derror[network.structinfo[offs+2]] = derror[network.structinfo[offs+2]]+dedf*dfdnet;
3110	}
3111	if( network.structinfo[offs+0]==0 )
3112	{
3113
3114	//
3115	// Adaptive summator
3116	//
3117	n1 = network.structinfo[offs+2];
3118	n2 = n1+network.structinfo[offs+1]-1;
3119	w1 = network.structinfo[offs+3];
3120	w2 = w1+network.structinfo[offs+1]-1;
3121	dedf = network.derror[i];
3122	dfdnet = 1.0;
3123	v = dedf*dfdnet;
3124	i1_ = (n1) - (w1);
3125	for(i_=w1; i_<=w2;i_++)
3126	{
3127	grad[i_] = v*neurons[i_+i1_];
3128	}
3129	i1_ = (w1) - (n1);
3130	for(i_=n1; i_<=n2;i_++)
3131	{
3132	derror[i_] = derror[i_] + v*weights[i_+i1_];
3133	}
3134	}
3135	if( network.structinfo[offs+0]<0 )
3136	{
3137	bflag = false;
3138	if( network.structinfo[offs+0]==-2 \| network.structinfo[offs+0]==-3 \| network.structinfo[offs+0]==-4 )
3139	{
3140
3141	//
3142	// Special neuron type, no back-propagation required
3143	//
3144	bflag = true;
3145	}
3146	System.Diagnostics.Debug.Assert(bflag, "MLPInternalCalculateGradient: unknown neuron type!");
3147	}
3148	}
3149	}
3150
3151
3152	/*************************************************************************
3153	Internal subroutine, chunked gradient
3154	*************************************************************************/
3155	private static void mlpchunkedgradient(ref multilayerperceptron network,
3156	ref double[,] xy,
3157	int cstart,
3158	int csize,
3159	ref double e,
3160	ref double[] grad,
3161	bool naturalerrorfunc)
3162	{
3163	int i = 0;
3164	int j = 0;
3165	int k = 0;
3166	int kl = 0;
3167	int n1 = 0;
3168	int n2 = 0;
3169	int w1 = 0;
3170	int w2 = 0;
3171	int c1 = 0;
3172	int c2 = 0;
3173	int ntotal = 0;
3174	int nin = 0;
3175	int nout = 0;
3176	int offs = 0;
3177	double dedf = 0;
3178	double dfdnet = 0;
3179	double f = 0;
3180	double df = 0;
3181	double d2f = 0;
3182	double v = 0;
3183	double s = 0;
3184	double fown = 0;
3185	double deown = 0;
3186	double net = 0;
3187	double lnnet = 0;
3188	double mx = 0;
3189	bool bflag = new bool();
3190	int istart = 0;
3191	int ineurons = 0;
3192	int idfdnet = 0;
3193	int iderror = 0;
3194	int izeros = 0;
3195	int i_ = 0;
3196	int i1_ = 0;
3197
3198
3199	//
3200	// Read network geometry, prepare data
3201	//
3202	nin = network.structinfo[1];
3203	nout = network.structinfo[2];
3204	ntotal = network.structinfo[3];
3205	istart = network.structinfo[5];
3206	c1 = cstart;
3207	c2 = cstart+csize-1;
3208	ineurons = 0;
3209	idfdnet = ntotal;
3210	iderror = 2*ntotal;
3211	izeros = 3*ntotal;
3212	for(j=0; j<=csize-1; j++)
3213	{
3214	network.chunks[izeros,j] = 0;
3215	}
3216
3217	//
3218	// Forward pass:
3219	// 1. Load inputs from XY to Chunks[0:NIn-1,0:CSize-1]
3220	// 2. Forward pass
3221	//
3222	for(i=0; i<=nin-1; i++)
3223	{
3224	for(j=0; j<=csize-1; j++)
3225	{
3226	if( (double)(network.columnsigmas[i])!=(double)(0) )
3227	{
3228	network.chunks[i,j] = (xy[c1+j,i]-network.columnmeans[i])/network.columnsigmas[i];
3229	}
3230	else
3231	{
3232	network.chunks[i,j] = xy[c1+j,i]-network.columnmeans[i];
3233	}
3234	}
3235	}
3236	for(i=0; i<=ntotal-1; i++)
3237	{
3238	offs = istart+i*nfieldwidth;
3239	if( network.structinfo[offs+0]>0 )
3240	{
3241
3242	//
3243	// Activation function:
3244	// * calculate F vector, F(i) = F(NET(i))
3245	//
3246	n1 = network.structinfo[offs+2];
3247	for(i_=0; i_<=csize-1;i_++)
3248	{
3249	network.chunks[i,i_] = network.chunks[n1,i_];
3250	}
3251	for(j=0; j<=csize-1; j++)
3252	{
3253	mlpactivationfunction(network.chunks[i,j], network.structinfo[offs+0], ref f, ref df, ref d2f);
3254	network.chunks[i,j] = f;
3255	network.chunks[idfdnet+i,j] = df;
3256	}
3257	}
3258	if( network.structinfo[offs+0]==0 )
3259	{
3260
3261	//
3262	// Adaptive summator:
3263	// * calculate NET vector, NET(i) = SUM(W(j,i)*Neurons(j),j=N1..N2)
3264	//
3265	n1 = network.structinfo[offs+2];
3266	n2 = n1+network.structinfo[offs+1]-1;
3267	w1 = network.structinfo[offs+3];
3268	w2 = w1+network.structinfo[offs+1]-1;
3269	for(i_=0; i_<=csize-1;i_++)
3270	{
3271	network.chunks[i,i_] = network.chunks[izeros,i_];
3272	}
3273	for(j=n1; j<=n2; j++)
3274	{
3275	v = network.weights[w1+j-n1];
3276	for(i_=0; i_<=csize-1;i_++)
3277	{
3278	network.chunks[i,i_] = network.chunks[i,i_] + v*network.chunks[j,i_];
3279	}
3280	}
3281	}
3282	if( network.structinfo[offs+0]<0 )
3283	{
3284	bflag = false;
3285	if( network.structinfo[offs+0]==-2 )
3286	{
3287
3288	//
3289	// input neuron, left unchanged
3290	//
3291	bflag = true;
3292	}
3293	if( network.structinfo[offs+0]==-3 )
3294	{
3295
3296	//
3297	// "-1" neuron
3298	//
3299	for(k=0; k<=csize-1; k++)
3300	{
3301	network.chunks[i,k] = -1;
3302	}
3303	bflag = true;
3304	}
3305	if( network.structinfo[offs+0]==-4 )
3306	{
3307
3308	//
3309	// "0" neuron
3310	//
3311	for(k=0; k<=csize-1; k++)
3312	{
3313	network.chunks[i,k] = 0;
3314	}
3315	bflag = true;
3316	}
3317	System.Diagnostics.Debug.Assert(bflag, "MLPChunkedGradient: internal error - unknown neuron type!");
3318	}
3319	}
3320
3321	//
3322	// Post-processing, error, dError/dOut
3323	//
3324	for(i=0; i<=ntotal-1; i++)
3325	{
3326	for(i_=0; i_<=csize-1;i_++)
3327	{
3328	network.chunks[iderror+i,i_] = network.chunks[izeros,i_];
3329	}
3330	}
3331	System.Diagnostics.Debug.Assert(network.structinfo[6]==0 \| network.structinfo[6]==1, "MLPChunkedGradient: unknown normalization type!");
3332	if( network.structinfo[6]==1 )
3333	{
3334
3335	//
3336	// Softmax output, classification network.
3337	//
3338	// For each K = 0..CSize-1 do:
3339	// 1. place exp(outputs[k]) to NWBuf[0:NOut-1]
3340	// 2. place sum(exp(..)) to NET
3341	// 3. calculate dError/dOut and place it to the second block of Chunks
3342	//
3343	for(k=0; k<=csize-1; k++)
3344	{
3345
3346	//
3347	// Normalize
3348	//
3349	mx = network.chunks[ntotal-nout,k];
3350	for(i=1; i<=nout-1; i++)
3351	{
3352	mx = Math.Max(mx, network.chunks[ntotal-nout+i,k]);
3353	}
3354	net = 0;
3355	for(i=0; i<=nout-1; i++)
3356	{
3357	network.nwbuf[i] = Math.Exp(network.chunks[ntotal-nout+i,k]-mx);
3358	net = net+network.nwbuf[i];
3359	}
3360
3361	//
3362	// Calculate error function and dError/dOut
3363	//
3364	if( naturalerrorfunc )
3365	{
3366
3367	//
3368	// Natural error func.
3369	//
3370	//
3371	s = 1;
3372	lnnet = Math.Log(net);
3373	kl = (int)Math.Round(xy[cstart+k,nin]);
3374	for(i=0; i<=nout-1; i++)
3375	{
3376	if( i==kl )
3377	{
3378	v = 1;
3379	}
3380	else
3381	{
3382	v = 0;
3383	}
3384	network.chunks[iderror+ntotal-nout+i,k] = s*network.nwbuf[i]/net-v;
3385	e = e+safecrossentropy(v, network.nwbuf[i]/net);
3386	}
3387	}
3388	else
3389	{
3390
3391	//
3392	// Least squares error func
3393	// Error, dError/dOut(normalized)
3394	//
3395	kl = (int)Math.Round(xy[cstart+k,nin]);
3396	for(i=0; i<=nout-1; i++)
3397	{
3398	if( i==kl )
3399	{
3400	v = network.nwbuf[i]/net-1;
3401	}
3402	else
3403	{
3404	v = network.nwbuf[i]/net;
3405	}
3406	network.nwbuf[nout+i] = v;
3407	e = e+AP.Math.Sqr(v)/2;
3408	}
3409
3410	//
3411	// From dError/dOut(normalized) to dError/dOut(non-normalized)
3412	//
3413	i1_ = (0)-(nout);
3414	v = 0.0;
3415	for(i_=nout; i_<=2*nout-1;i_++)
3416	{
3417	v += network.nwbuf[i_]*network.nwbuf[i_+i1_];
3418	}
3419	for(i=0; i<=nout-1; i++)
3420	{
3421	fown = network.nwbuf[i];
3422	deown = network.nwbuf[nout+i];
3423	network.chunks[iderror+ntotal-nout+i,k] = (-v+deownfown+deown(net-fown))*fown/AP.Math.Sqr(net);
3424	}
3425	}
3426	}
3427	}
3428	else
3429	{
3430
3431	//
3432	// Normal output, regression network
3433	//
3434	// For each K = 0..CSize-1 do:
3435	// 1. calculate dError/dOut and place it to the second block of Chunks
3436	//
3437	for(i=0; i<=nout-1; i++)
3438	{
3439	for(j=0; j<=csize-1; j++)
3440	{
3441	v = network.chunks[ntotal-nout+i,j]*network.columnsigmas[nin+i]+network.columnmeans[nin+i]-xy[cstart+j,nin+i];
3442	network.chunks[iderror+ntotal-nout+i,j] = v*network.columnsigmas[nin+i];
3443	e = e+AP.Math.Sqr(v)/2;
3444	}
3445	}
3446	}
3447
3448	//
3449	// Backpropagation
3450	//
3451	for(i=ntotal-1; i>=0; i--)
3452	{
3453
3454	//
3455	// Extract info
3456	//
3457	offs = istart+i*nfieldwidth;
3458	if( network.structinfo[offs+0]>0 )
3459	{
3460
3461	//
3462	// Activation function
3463	//
3464	n1 = network.structinfo[offs+2];
3465	for(k=0; k<=csize-1; k++)
3466	{
3467	network.chunks[iderror+i,k] = network.chunks[iderror+i,k]*network.chunks[idfdnet+i,k];
3468	}
3469	for(i_=0; i_<=csize-1;i_++)
3470	{
3471	network.chunks[iderror+n1,i_] = network.chunks[iderror+n1,i_] + network.chunks[iderror+i,i_];
3472	}
3473	}
3474	if( network.structinfo[offs+0]==0 )
3475	{
3476
3477	//
3478	// "Normal" activation function
3479	//
3480	n1 = network.structinfo[offs+2];
3481	n2 = n1+network.structinfo[offs+1]-1;
3482	w1 = network.structinfo[offs+3];
3483	w2 = w1+network.structinfo[offs+1]-1;
3484	for(j=w1; j<=w2; j++)
3485	{
3486	v = 0.0;
3487	for(i_=0; i_<=csize-1;i_++)
3488	{
3489	v += network.chunks[n1+j-w1,i_]*network.chunks[iderror+i,i_];
3490	}
3491	grad[j] = grad[j]+v;
3492	}
3493	for(j=n1; j<=n2; j++)
3494	{
3495	v = network.weights[w1+j-n1];
3496	for(i_=0; i_<=csize-1;i_++)
3497	{
3498	network.chunks[iderror+j,i_] = network.chunks[iderror+j,i_] + v*network.chunks[iderror+i,i_];
3499	}
3500	}
3501	}
3502	if( network.structinfo[offs+0]<0 )
3503	{
3504	bflag = false;
3505	if( network.structinfo[offs+0]==-2 \| network.structinfo[offs+0]==-3 \| network.structinfo[offs+0]==-4 )
3506	{
3507
3508	//
3509	// Special neuron type, no back-propagation required
3510	//
3511	bflag = true;
3512	}
3513	System.Diagnostics.Debug.Assert(bflag, "MLPInternalCalculateGradient: unknown neuron type!");
3514	}
3515	}
3516	}
3517
3518
3519	/*************************************************************************
3520	Returns T*Ln(T/Z), guarded against overflow/underflow.
3521	Internal subroutine.
3522	*************************************************************************/
3523	private static double safecrossentropy(double t,
3524	double z)
3525	{
3526	double result = 0;
3527	double r = 0;
3528
3529	if( (double)(t)==(double)(0) )
3530	{
3531	result = 0;
3532	}
3533	else
3534	{
3535	if( (double)(Math.Abs(z))>(double)(1) )
3536	{
3537
3538	//
3539	// Shouldn't be the case with softmax,
3540	// but we just want to be sure.
3541	//
3542	if( (double)(t/z)==(double)(0) )
3543	{
3544	r = AP.Math.MinRealNumber;
3545	}
3546	else
3547	{
3548	r = t/z;
3549	}
3550	}
3551	else
3552	{
3553
3554	//
3555	// Normal case
3556	//
3557	if( (double)(z)==(double)(0) \| (double)(Math.Abs(t))>=(double)(AP.Math.MaxRealNumber*Math.Abs(z)) )
3558	{
3559	r = AP.Math.MaxRealNumber;
3560	}
3561	else
3562	{
3563	r = t/z;
3564	}
3565	}
3566	result = t*Math.Log(r);
3567	}
3568	return result;
3569	}
3570	}
3571	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats:

Update cookies preferences