Context Navigation

source: trunk/sources/HeuristicLab.ExtLibs/HeuristicLab.ALGLIB/2.5.0/ALGLIB-2.5.0/mlpbase.cs @ 3932

Visit:

Last change on this file since 3932 was 3839, checked in by mkommend, 15 years ago
implemented first version of LR (ticket #1012)
File size: 126.4 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 k = 0;
1327	int nin = 0;
1328	int nout = 0;
1329	int wcount = 0;
1330	double e = 0;
1331	int i_ = 0;
1332	int i1_ = 0;
1333
1334	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1335	result = 0;
1336	for(i=0; i<=ssize-1; i++)
1337	{
1338
1339	//
1340	// Process vector
1341	//
1342	for(i_=0; i_<=nin-1;i_++)
1343	{
1344	network.x[i_] = xy[i,i_];
1345	}
1346	mlpprocess(ref network, ref network.x, ref network.y);
1347
1348	//
1349	// Update error function
1350	//
1351	if( network.structinfo[6]==0 )
1352	{
1353
1354	//
1355	// Least squares error function
1356	//
1357	i1_ = (nin) - (0);
1358	for(i_=0; i_<=nout-1;i_++)
1359	{
1360	network.y[i_] = network.y[i_] - xy[i,i_+i1_];
1361	}
1362	e = 0.0;
1363	for(i_=0; i_<=nout-1;i_++)
1364	{
1365	e += network.y[i_]*network.y[i_];
1366	}
1367	result = result+e/2;
1368	}
1369	else
1370	{
1371
1372	//
1373	// Cross-entropy error function
1374	//
1375	k = (int)Math.Round(xy[i,nin]);
1376	if( k>=0 & k<nout )
1377	{
1378	result = result+safecrossentropy(1, network.y[k]);
1379	}
1380	}
1381	}
1382	return result;
1383	}
1384
1385
1386	/*************************************************************************
1387	Classification error
1388
1389	-- ALGLIB --
1390	Copyright 04.11.2007 by Bochkanov Sergey
1391	*************************************************************************/
1392	public static int mlpclserror(ref multilayerperceptron network,
1393	ref double[,] xy,
1394	int ssize)
1395	{
1396	int result = 0;
1397	int i = 0;
1398	int j = 0;
1399	int nin = 0;
1400	int nout = 0;
1401	int wcount = 0;
1402	double[] workx = new double[0];
1403	double[] worky = new double[0];
1404	int nn = 0;
1405	int ns = 0;
1406	int nmax = 0;
1407	int i_ = 0;
1408
1409	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1410	workx = new double[nin-1+1];
1411	worky = new double[nout-1+1];
1412	result = 0;
1413	for(i=0; i<=ssize-1; i++)
1414	{
1415
1416	//
1417	// Process
1418	//
1419	for(i_=0; i_<=nin-1;i_++)
1420	{
1421	workx[i_] = xy[i,i_];
1422	}
1423	mlpprocess(ref network, ref workx, ref worky);
1424
1425	//
1426	// Network version of the answer
1427	//
1428	nmax = 0;
1429	for(j=0; j<=nout-1; j++)
1430	{
1431	if( (double)(worky[j])>(double)(worky[nmax]) )
1432	{
1433	nmax = j;
1434	}
1435	}
1436	nn = nmax;
1437
1438	//
1439	// Right answer
1440	//
1441	if( mlpissoftmax(ref network) )
1442	{
1443	ns = (int)Math.Round(xy[i,nin]);
1444	}
1445	else
1446	{
1447	nmax = 0;
1448	for(j=0; j<=nout-1; j++)
1449	{
1450	if( (double)(xy[i,nin+j])>(double)(xy[i,nin+nmax]) )
1451	{
1452	nmax = j;
1453	}
1454	}
1455	ns = nmax;
1456	}
1457
1458	//
1459	// compare
1460	//
1461	if( nn!=ns )
1462	{
1463	result = result+1;
1464	}
1465	}
1466	return result;
1467	}
1468
1469
1470	/*************************************************************************
1471	Relative classification error on the test set
1472
1473	INPUT PARAMETERS:
1474	Network - network
1475	XY - test set
1476	NPoints - test set size
1477
1478	RESULT:
1479	percent of incorrectly classified cases. Works both for
1480	classifier networks and general purpose networks used as
1481	classifiers.
1482
1483	-- ALGLIB --
1484	Copyright 25.12.2008 by Bochkanov Sergey
1485	*************************************************************************/
1486	public static double mlprelclserror(ref multilayerperceptron network,
1487	ref double[,] xy,
1488	int npoints)
1489	{
1490	double result = 0;
1491
1492	result = (double)(mlpclserror(ref network, ref xy, npoints))/(double)(npoints);
1493	return result;
1494	}
1495
1496
1497	/*************************************************************************
1498	Average cross-entropy (in bits per element) on the test set
1499
1500	INPUT PARAMETERS:
1501	Network - neural network
1502	XY - test set
1503	NPoints - test set size
1504
1505	RESULT:
1506	CrossEntropy/(NPoints*LN(2)).
1507	Zero if network solves regression task.
1508
1509	-- ALGLIB --
1510	Copyright 08.01.2009 by Bochkanov Sergey
1511	*************************************************************************/
1512	public static double mlpavgce(ref multilayerperceptron network,
1513	ref double[,] xy,
1514	int npoints)
1515	{
1516	double result = 0;
1517	int nin = 0;
1518	int nout = 0;
1519	int wcount = 0;
1520
1521	if( mlpissoftmax(ref network) )
1522	{
1523	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1524	result = mlperrorn(ref network, ref xy, npoints)/(npoints*Math.Log(2));
1525	}
1526	else
1527	{
1528	result = 0;
1529	}
1530	return result;
1531	}
1532
1533
1534	/*************************************************************************
1535	RMS error on the test set
1536
1537	INPUT PARAMETERS:
1538	Network - neural network
1539	XY - test set
1540	NPoints - test set size
1541
1542	RESULT:
1543	root mean square error.
1544	Its meaning for regression task is obvious. As for
1545	classification task, RMS error means error when estimating posterior
1546	probabilities.
1547
1548	-- ALGLIB --
1549	Copyright 04.11.2007 by Bochkanov Sergey
1550	*************************************************************************/
1551	public static double mlprmserror(ref multilayerperceptron network,
1552	ref double[,] xy,
1553	int npoints)
1554	{
1555	double result = 0;
1556	int nin = 0;
1557	int nout = 0;
1558	int wcount = 0;
1559
1560	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1561	result = Math.Sqrt(2mlperror(ref network, ref xy, npoints)/(npointsnout));
1562	return result;
1563	}
1564
1565
1566	/*************************************************************************
1567	Average error on the test set
1568
1569	INPUT PARAMETERS:
1570	Network - neural network
1571	XY - test set
1572	NPoints - test set size
1573
1574	RESULT:
1575	Its meaning for regression task is obvious. As for
1576	classification task, it means average error when estimating posterior
1577	probabilities.
1578
1579	-- ALGLIB --
1580	Copyright 11.03.2008 by Bochkanov Sergey
1581	*************************************************************************/
1582	public static double mlpavgerror(ref multilayerperceptron network,
1583	ref double[,] xy,
1584	int npoints)
1585	{
1586	double result = 0;
1587	int i = 0;
1588	int j = 0;
1589	int k = 0;
1590	int nin = 0;
1591	int nout = 0;
1592	int wcount = 0;
1593	int i_ = 0;
1594
1595	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1596	result = 0;
1597	for(i=0; i<=npoints-1; i++)
1598	{
1599	for(i_=0; i_<=nin-1;i_++)
1600	{
1601	network.x[i_] = xy[i,i_];
1602	}
1603	mlpprocess(ref network, ref network.x, ref network.y);
1604	if( mlpissoftmax(ref network) )
1605	{
1606
1607	//
1608	// class labels
1609	//
1610	k = (int)Math.Round(xy[i,nin]);
1611	for(j=0; j<=nout-1; j++)
1612	{
1613	if( j==k )
1614	{
1615	result = result+Math.Abs(1-network.y[j]);
1616	}
1617	else
1618	{
1619	result = result+Math.Abs(network.y[j]);
1620	}
1621	}
1622	}
1623	else
1624	{
1625
1626	//
1627	// real outputs
1628	//
1629	for(j=0; j<=nout-1; j++)
1630	{
1631	result = result+Math.Abs(xy[i,nin+j]-network.y[j]);
1632	}
1633	}
1634	}
1635	result = result/(npoints*nout);
1636	return result;
1637	}
1638
1639
1640	/*************************************************************************
1641	Average relative error on the test set
1642
1643	INPUT PARAMETERS:
1644	Network - neural network
1645	XY - test set
1646	NPoints - test set size
1647
1648	RESULT:
1649	Its meaning for regression task is obvious. As for
1650	classification task, it means average relative error when estimating
1651	posterior probability of belonging to the correct class.
1652
1653	-- ALGLIB --
1654	Copyright 11.03.2008 by Bochkanov Sergey
1655	*************************************************************************/
1656	public static double mlpavgrelerror(ref multilayerperceptron network,
1657	ref double[,] xy,
1658	int npoints)
1659	{
1660	double result = 0;
1661	int i = 0;
1662	int j = 0;
1663	int k = 0;
1664	int lk = 0;
1665	int nin = 0;
1666	int nout = 0;
1667	int wcount = 0;
1668	int i_ = 0;
1669
1670	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1671	result = 0;
1672	k = 0;
1673	for(i=0; i<=npoints-1; i++)
1674	{
1675	for(i_=0; i_<=nin-1;i_++)
1676	{
1677	network.x[i_] = xy[i,i_];
1678	}
1679	mlpprocess(ref network, ref network.x, ref network.y);
1680	if( mlpissoftmax(ref network) )
1681	{
1682
1683	//
1684	// class labels
1685	//
1686	lk = (int)Math.Round(xy[i,nin]);
1687	for(j=0; j<=nout-1; j++)
1688	{
1689	if( j==lk )
1690	{
1691	result = result+Math.Abs(1-network.y[j]);
1692	k = k+1;
1693	}
1694	}
1695	}
1696	else
1697	{
1698
1699	//
1700	// real outputs
1701	//
1702	for(j=0; j<=nout-1; j++)
1703	{
1704	if( (double)(xy[i,nin+j])!=(double)(0) )
1705	{
1706	result = result+Math.Abs(xy[i,nin+j]-network.y[j])/Math.Abs(xy[i,nin+j]);
1707	k = k+1;
1708	}
1709	}
1710	}
1711	}
1712	if( k!=0 )
1713	{
1714	result = result/k;
1715	}
1716	return result;
1717	}
1718
1719
1720	/*************************************************************************
1721	Gradient calculation. Internal subroutine.
1722
1723	-- ALGLIB --
1724	Copyright 04.11.2007 by Bochkanov Sergey
1725	*************************************************************************/
1726	public static void mlpgrad(ref multilayerperceptron network,
1727	ref double[] x,
1728	ref double[] desiredy,
1729	ref double e,
1730	ref double[] grad)
1731	{
1732	int i = 0;
1733	int nout = 0;
1734	int ntotal = 0;
1735
1736
1737	//
1738	// Prepare dError/dOut, internal structures
1739	//
1740	mlpprocess(ref network, ref x, ref network.y);
1741	nout = network.structinfo[2];
1742	ntotal = network.structinfo[3];
1743	e = 0;
1744	for(i=0; i<=ntotal-1; i++)
1745	{
1746	network.derror[i] = 0;
1747	}
1748	for(i=0; i<=nout-1; i++)
1749	{
1750	network.derror[ntotal-nout+i] = network.y[i]-desiredy[i];
1751	e = e+AP.Math.Sqr(network.y[i]-desiredy[i])/2;
1752	}
1753
1754	//
1755	// gradient
1756	//
1757	mlpinternalcalculategradient(ref network, ref network.neurons, ref network.weights, ref network.derror, ref grad, false);
1758	}
1759
1760
1761	/*************************************************************************
1762	Gradient calculation (natural error function). Internal subroutine.
1763
1764	-- ALGLIB --
1765	Copyright 04.11.2007 by Bochkanov Sergey
1766	*************************************************************************/
1767	public static void mlpgradn(ref multilayerperceptron network,
1768	ref double[] x,
1769	ref double[] desiredy,
1770	ref double e,
1771	ref double[] grad)
1772	{
1773	double s = 0;
1774	int i = 0;
1775	int nout = 0;
1776	int ntotal = 0;
1777
1778
1779	//
1780	// Prepare dError/dOut, internal structures
1781	//
1782	mlpprocess(ref network, ref x, ref network.y);
1783	nout = network.structinfo[2];
1784	ntotal = network.structinfo[3];
1785	for(i=0; i<=ntotal-1; i++)
1786	{
1787	network.derror[i] = 0;
1788	}
1789	e = 0;
1790	if( network.structinfo[6]==0 )
1791	{
1792
1793	//
1794	// Regression network, least squares
1795	//
1796	for(i=0; i<=nout-1; i++)
1797	{
1798	network.derror[ntotal-nout+i] = network.y[i]-desiredy[i];
1799	e = e+AP.Math.Sqr(network.y[i]-desiredy[i])/2;
1800	}
1801	}
1802	else
1803	{
1804
1805	//
1806	// Classification network, cross-entropy
1807	//
1808	s = 0;
1809	for(i=0; i<=nout-1; i++)
1810	{
1811	s = s+desiredy[i];
1812	}
1813	for(i=0; i<=nout-1; i++)
1814	{
1815	network.derror[ntotal-nout+i] = s*network.y[i]-desiredy[i];
1816	e = e+safecrossentropy(desiredy[i], network.y[i]);
1817	}
1818	}
1819
1820	//
1821	// gradient
1822	//
1823	mlpinternalcalculategradient(ref network, ref network.neurons, ref network.weights, ref network.derror, ref grad, true);
1824	}
1825
1826
1827	/*************************************************************************
1828	Batch gradient calculation. Internal subroutine.
1829
1830	-- ALGLIB --
1831	Copyright 04.11.2007 by Bochkanov Sergey
1832	*************************************************************************/
1833	public static void mlpgradbatch(ref multilayerperceptron network,
1834	ref double[,] xy,
1835	int ssize,
1836	ref double e,
1837	ref double[] grad)
1838	{
1839	int i = 0;
1840	int nin = 0;
1841	int nout = 0;
1842	int wcount = 0;
1843
1844	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1845	for(i=0; i<=wcount-1; i++)
1846	{
1847	grad[i] = 0;
1848	}
1849	e = 0;
1850	i = 0;
1851	while( i<=ssize-1 )
1852	{
1853	mlpchunkedgradient(ref network, ref xy, i, Math.Min(ssize, i+chunksize)-i, ref e, ref grad, false);
1854	i = i+chunksize;
1855	}
1856	}
1857
1858
1859	/*************************************************************************
1860	Batch gradient calculation (natural error function). Internal subroutine.
1861
1862	-- ALGLIB --
1863	Copyright 04.11.2007 by Bochkanov Sergey
1864	*************************************************************************/
1865	public static void mlpgradnbatch(ref multilayerperceptron network,
1866	ref double[,] xy,
1867	int ssize,
1868	ref double e,
1869	ref double[] grad)
1870	{
1871	int i = 0;
1872	int nin = 0;
1873	int nout = 0;
1874	int wcount = 0;
1875
1876	mlpproperties(ref network, ref nin, ref nout, ref wcount);
1877	for(i=0; i<=wcount-1; i++)
1878	{
1879	grad[i] = 0;
1880	}
1881	e = 0;
1882	i = 0;
1883	while( i<=ssize-1 )
1884	{
1885	mlpchunkedgradient(ref network, ref xy, i, Math.Min(ssize, i+chunksize)-i, ref e, ref grad, true);
1886	i = i+chunksize;
1887	}
1888	}
1889
1890
1891	/*************************************************************************
1892	Batch Hessian calculation (natural error function) using R-algorithm.
1893	Internal subroutine.
1894
1895	-- ALGLIB --
1896	Copyright 26.01.2008 by Bochkanov Sergey.
1897
1898	Hessian calculation based on R-algorithm described in
1899	"Fast Exact Multiplication by the Hessian",
1900	B. A. Pearlmutter,
1901	Neural Computation, 1994.
1902	*************************************************************************/
1903	public static void mlphessiannbatch(ref multilayerperceptron network,
1904	ref double[,] xy,
1905	int ssize,
1906	ref double e,
1907	ref double[] grad,
1908	ref double[,] h)
1909	{
1910	mlphessianbatchinternal(ref network, ref xy, ssize, true, ref e, ref grad, ref h);
1911	}
1912
1913
1914	/*************************************************************************
1915	Batch Hessian calculation using R-algorithm.
1916	Internal subroutine.
1917
1918	-- ALGLIB --
1919	Copyright 26.01.2008 by Bochkanov Sergey.
1920
1921	Hessian calculation based on R-algorithm described in
1922	"Fast Exact Multiplication by the Hessian",
1923	B. A. Pearlmutter,
1924	Neural Computation, 1994.
1925	*************************************************************************/
1926	public static void mlphessianbatch(ref multilayerperceptron network,
1927	ref double[,] xy,
1928	int ssize,
1929	ref double e,
1930	ref double[] grad,
1931	ref double[,] h)
1932	{
1933	mlphessianbatchinternal(ref network, ref xy, ssize, false, ref e, ref grad, ref h);
1934	}
1935
1936
1937	/*************************************************************************
1938	Internal subroutine, shouldn't be called by user.
1939	*************************************************************************/
1940	public static void mlpinternalprocessvector(ref int[] structinfo,
1941	ref double[] weights,
1942	ref double[] columnmeans,
1943	ref double[] columnsigmas,
1944	ref double[] neurons,
1945	ref double[] dfdnet,
1946	ref double[] x,
1947	ref double[] y)
1948	{
1949	int i = 0;
1950	int n1 = 0;
1951	int n2 = 0;
1952	int w1 = 0;
1953	int w2 = 0;
1954	int ntotal = 0;
1955	int nin = 0;
1956	int nout = 0;
1957	int istart = 0;
1958	int offs = 0;
1959	double net = 0;
1960	double f = 0;
1961	double df = 0;
1962	double d2f = 0;
1963	double mx = 0;
1964	bool perr = new bool();
1965	int i_ = 0;
1966	int i1_ = 0;
1967
1968
1969	//
1970	// Read network geometry
1971	//
1972	nin = structinfo[1];
1973	nout = structinfo[2];
1974	ntotal = structinfo[3];
1975	istart = structinfo[5];
1976
1977	//
1978	// Inputs standartisation and putting in the network
1979	//
1980	for(i=0; i<=nin-1; i++)
1981	{
1982	if( (double)(columnsigmas[i])!=(double)(0) )
1983	{
1984	neurons[i] = (x[i]-columnmeans[i])/columnsigmas[i];
1985	}
1986	else
1987	{
1988	neurons[i] = x[i]-columnmeans[i];
1989	}
1990	}
1991
1992	//
1993	// Process network
1994	//
1995	for(i=0; i<=ntotal-1; i++)
1996	{
1997	offs = istart+i*nfieldwidth;
1998	if( structinfo[offs+0]>0 )
1999	{
2000
2001	//
2002	// Activation function
2003	//
2004	mlpactivationfunction(neurons[structinfo[offs+2]], structinfo[offs+0], ref f, ref df, ref d2f);
2005	neurons[i] = f;
2006	dfdnet[i] = df;
2007	}
2008	if( structinfo[offs+0]==0 )
2009	{
2010
2011	//
2012	// Adaptive summator
2013	//
2014	n1 = structinfo[offs+2];
2015	n2 = n1+structinfo[offs+1]-1;
2016	w1 = structinfo[offs+3];
2017	w2 = w1+structinfo[offs+1]-1;
2018	i1_ = (n1)-(w1);
2019	net = 0.0;
2020	for(i_=w1; i_<=w2;i_++)
2021	{
2022	net += weights[i_]*neurons[i_+i1_];
2023	}
2024	neurons[i] = net;
2025	dfdnet[i] = 1.0;
2026	}
2027	if( structinfo[offs+0]<0 )
2028	{
2029	perr = true;
2030	if( structinfo[offs+0]==-2 )
2031	{
2032
2033	//
2034	// input neuron, left unchanged
2035	//
2036	perr = false;
2037	}
2038	if( structinfo[offs+0]==-3 )
2039	{
2040
2041	//
2042	// "-1" neuron
2043	//
2044	neurons[i] = -1;
2045	perr = false;
2046	}
2047	if( structinfo[offs+0]==-4 )
2048	{
2049
2050	//
2051	// "0" neuron
2052	//
2053	neurons[i] = 0;
2054	perr = false;
2055	}
2056	System.Diagnostics.Debug.Assert(!perr, "MLPInternalProcessVector: internal error - unknown neuron type!");
2057	}
2058	}
2059
2060	//
2061	// Extract result
2062	//
2063	i1_ = (ntotal-nout) - (0);
2064	for(i_=0; i_<=nout-1;i_++)
2065	{
2066	y[i_] = neurons[i_+i1_];
2067	}
2068
2069	//
2070	// Softmax post-processing or standardisation if needed
2071	//
2072	System.Diagnostics.Debug.Assert(structinfo[6]==0 \| structinfo[6]==1, "MLPInternalProcessVector: unknown normalization type!");
2073	if( structinfo[6]==1 )
2074	{
2075
2076	//
2077	// Softmax
2078	//
2079	mx = y[0];
2080	for(i=1; i<=nout-1; i++)
2081	{
2082	mx = Math.Max(mx, y[i]);
2083	}
2084	net = 0;
2085	for(i=0; i<=nout-1; i++)
2086	{
2087	y[i] = Math.Exp(y[i]-mx);
2088	net = net+y[i];
2089	}
2090	for(i=0; i<=nout-1; i++)
2091	{
2092	y[i] = y[i]/net;
2093	}
2094	}
2095	else
2096	{
2097
2098	//
2099	// Standardisation
2100	//
2101	for(i=0; i<=nout-1; i++)
2102	{
2103	y[i] = y[i]*columnsigmas[nin+i]+columnmeans[nin+i];
2104	}
2105	}
2106	}
2107
2108
2109	/*************************************************************************
2110	Internal subroutine: adding new input layer to network
2111	*************************************************************************/
2112	private static void addinputlayer(int ncount,
2113	ref int[] lsizes,
2114	ref int[] ltypes,
2115	ref int[] lconnfirst,
2116	ref int[] lconnlast,
2117	ref int lastproc)
2118	{
2119	lsizes[0] = ncount;
2120	ltypes[0] = -2;
2121	lconnfirst[0] = 0;
2122	lconnlast[0] = 0;
2123	lastproc = 0;
2124	}
2125
2126
2127	/*************************************************************************
2128	Internal subroutine: adding new summator layer to network
2129	*************************************************************************/
2130	private static void addbiasedsummatorlayer(int ncount,
2131	ref int[] lsizes,
2132	ref int[] ltypes,
2133	ref int[] lconnfirst,
2134	ref int[] lconnlast,
2135	ref int lastproc)
2136	{
2137	lsizes[lastproc+1] = 1;
2138	ltypes[lastproc+1] = -3;
2139	lconnfirst[lastproc+1] = 0;
2140	lconnlast[lastproc+1] = 0;
2141	lsizes[lastproc+2] = ncount;
2142	ltypes[lastproc+2] = 0;
2143	lconnfirst[lastproc+2] = lastproc;
2144	lconnlast[lastproc+2] = lastproc+1;
2145	lastproc = lastproc+2;
2146	}
2147
2148
2149	/*************************************************************************
2150	Internal subroutine: adding new summator layer to network
2151	*************************************************************************/
2152	private static void addactivationlayer(int functype,
2153	ref int[] lsizes,
2154	ref int[] ltypes,
2155	ref int[] lconnfirst,
2156	ref int[] lconnlast,
2157	ref int lastproc)
2158	{
2159	System.Diagnostics.Debug.Assert(functype>0, "AddActivationLayer: incorrect function type");
2160	lsizes[lastproc+1] = lsizes[lastproc];
2161	ltypes[lastproc+1] = functype;
2162	lconnfirst[lastproc+1] = lastproc;
2163	lconnlast[lastproc+1] = lastproc;
2164	lastproc = lastproc+1;
2165	}
2166
2167
2168	/*************************************************************************
2169	Internal subroutine: adding new zero layer to network
2170	*************************************************************************/
2171	private static void addzerolayer(ref int[] lsizes,
2172	ref int[] ltypes,
2173	ref int[] lconnfirst,
2174	ref int[] lconnlast,
2175	ref int lastproc)
2176	{
2177	lsizes[lastproc+1] = 1;
2178	ltypes[lastproc+1] = -4;
2179	lconnfirst[lastproc+1] = 0;
2180	lconnlast[lastproc+1] = 0;
2181	lastproc = lastproc+1;
2182	}
2183
2184
2185	/*************************************************************************
2186	Internal subroutine.
2187
2188	-- ALGLIB --
2189	Copyright 04.11.2007 by Bochkanov Sergey
2190	*************************************************************************/
2191	private static void mlpcreate(int nin,
2192	int nout,
2193	ref int[] lsizes,
2194	ref int[] ltypes,
2195	ref int[] lconnfirst,
2196	ref int[] lconnlast,
2197	int layerscount,
2198	bool isclsnet,
2199	ref multilayerperceptron network)
2200	{
2201	int i = 0;
2202	int j = 0;
2203	int ssize = 0;
2204	int ntotal = 0;
2205	int wcount = 0;
2206	int offs = 0;
2207	int nprocessed = 0;
2208	int wallocated = 0;
2209	int[] localtemp = new int[0];
2210	int[] lnfirst = new int[0];
2211	int[] lnsyn = new int[0];
2212
2213
2214	//
2215	// Check
2216	//
2217	System.Diagnostics.Debug.Assert(layerscount>0, "MLPCreate: wrong parameters!");
2218	System.Diagnostics.Debug.Assert(ltypes[0]==-2, "MLPCreate: wrong LTypes[0] (must be -2)!");
2219	for(i=0; i<=layerscount-1; i++)
2220	{
2221	System.Diagnostics.Debug.Assert(lsizes[i]>0, "MLPCreate: wrong LSizes!");
2222	System.Diagnostics.Debug.Assert(lconnfirst[i]>=0 & (lconnfirst[i]<i \| i==0), "MLPCreate: wrong LConnFirst!");
2223	System.Diagnostics.Debug.Assert(lconnlast[i]>=lconnfirst[i] & (lconnlast[i]<i \| i==0), "MLPCreate: wrong LConnLast!");
2224	}
2225
2226	//
2227	// Build network geometry
2228	//
2229	lnfirst = new int[layerscount-1+1];
2230	lnsyn = new int[layerscount-1+1];
2231	ntotal = 0;
2232	wcount = 0;
2233	for(i=0; i<=layerscount-1; i++)
2234	{
2235
2236	//
2237	// Analyze connections.
2238	// This code must throw an assertion in case of unknown LTypes[I]
2239	//
2240	lnsyn[i] = -1;
2241	if( ltypes[i]>=0 )
2242	{
2243	lnsyn[i] = 0;
2244	for(j=lconnfirst[i]; j<=lconnlast[i]; j++)
2245	{
2246	lnsyn[i] = lnsyn[i]+lsizes[j];
2247	}
2248	}
2249	else
2250	{
2251	if( ltypes[i]==-2 \| ltypes[i]==-3 \| ltypes[i]==-4 )
2252	{
2253	lnsyn[i] = 0;
2254	}
2255	}
2256	System.Diagnostics.Debug.Assert(lnsyn[i]>=0, "MLPCreate: internal error #0!");
2257
2258	//
2259	// Other info
2260	//
2261	lnfirst[i] = ntotal;
2262	ntotal = ntotal+lsizes[i];
2263	if( ltypes[i]==0 )
2264	{
2265	wcount = wcount+lnsyn[i]*lsizes[i];
2266	}
2267	}
2268	ssize = 7+ntotal*nfieldwidth;
2269
2270	//
2271	// Allocate
2272	//
2273	network.structinfo = new int[ssize-1+1];
2274	network.weights = new double[wcount-1+1];
2275	if( isclsnet )
2276	{
2277	network.columnmeans = new double[nin-1+1];
2278	network.columnsigmas = new double[nin-1+1];
2279	}
2280	else
2281	{
2282	network.columnmeans = new double[nin+nout-1+1];
2283	network.columnsigmas = new double[nin+nout-1+1];
2284	}
2285	network.neurons = new double[ntotal-1+1];
2286	network.chunks = new double[3*ntotal+1, chunksize-1+1];
2287	network.nwbuf = new double[Math.Max(wcount, 2*nout)-1+1];
2288	network.dfdnet = new double[ntotal-1+1];
2289	network.x = new double[nin-1+1];
2290	network.y = new double[nout-1+1];
2291	network.derror = new double[ntotal-1+1];
2292
2293	//
2294	// Fill structure: global info
2295	//
2296	network.structinfo[0] = ssize;
2297	network.structinfo[1] = nin;
2298	network.structinfo[2] = nout;
2299	network.structinfo[3] = ntotal;
2300	network.structinfo[4] = wcount;
2301	network.structinfo[5] = 7;
2302	if( isclsnet )
2303	{
2304	network.structinfo[6] = 1;
2305	}
2306	else
2307	{
2308	network.structinfo[6] = 0;
2309	}
2310
2311	//
2312	// Fill structure: neuron connections
2313	//
2314	nprocessed = 0;
2315	wallocated = 0;
2316	for(i=0; i<=layerscount-1; i++)
2317	{
2318	for(j=0; j<=lsizes[i]-1; j++)
2319	{
2320	offs = network.structinfo[5]+nprocessed*nfieldwidth;
2321	network.structinfo[offs+0] = ltypes[i];
2322	if( ltypes[i]==0 )
2323	{
2324
2325	//
2326	// Adaptive summator:
2327	// * connections with weights to previous neurons
2328	//
2329	network.structinfo[offs+1] = lnsyn[i];
2330	network.structinfo[offs+2] = lnfirst[lconnfirst[i]];
2331	network.structinfo[offs+3] = wallocated;
2332	wallocated = wallocated+lnsyn[i];
2333	nprocessed = nprocessed+1;
2334	}
2335	if( ltypes[i]>0 )
2336	{
2337
2338	//
2339	// Activation layer:
2340	// * each neuron connected to one (only one) of previous neurons.
2341	// * no weights
2342	//
2343	network.structinfo[offs+1] = 1;
2344	network.structinfo[offs+2] = lnfirst[lconnfirst[i]]+j;
2345	network.structinfo[offs+3] = -1;
2346	nprocessed = nprocessed+1;
2347	}
2348	if( ltypes[i]==-2 \| ltypes[i]==-3 \| ltypes[i]==-4 )
2349	{
2350	nprocessed = nprocessed+1;
2351	}
2352	}
2353	}
2354	System.Diagnostics.Debug.Assert(wallocated==wcount, "MLPCreate: internal error #1!");
2355	System.Diagnostics.Debug.Assert(nprocessed==ntotal, "MLPCreate: internal error #2!");
2356
2357	//
2358	// Fill weights by small random values
2359	// Initialize means and sigmas
2360	//
2361	for(i=0; i<=wcount-1; i++)
2362	{
2363	network.weights[i] = AP.Math.RandomReal()-0.5;
2364	}
2365	for(i=0; i<=nin-1; i++)
2366	{
2367	network.columnmeans[i] = 0;
2368	network.columnsigmas[i] = 1;
2369	}
2370	if( !isclsnet )
2371	{
2372	for(i=0; i<=nout-1; i++)
2373	{
2374	network.columnmeans[nin+i] = 0;
2375	network.columnsigmas[nin+i] = 1;
2376	}
2377	}
2378	}
2379
2380
2381	/*************************************************************************
2382	Internal subroutine
2383
2384	-- ALGLIB --
2385	Copyright 04.11.2007 by Bochkanov Sergey
2386	*************************************************************************/
2387	private static void mlpactivationfunction(double net,
2388	int k,
2389	ref double f,
2390	ref double df,
2391	ref double d2f)
2392	{
2393	double net2 = 0;
2394	double arg = 0;
2395	double root = 0;
2396	double r = 0;
2397
2398	f = 0;
2399	df = 0;
2400	if( k==1 )
2401	{
2402
2403	//
2404	// TanH activation function
2405	//
2406	if( (double)(Math.Abs(net))<(double)(100) )
2407	{
2408	f = Math.Tanh(net);
2409	}
2410	else
2411	{
2412	f = Math.Sign(net);
2413	}
2414	df = 1-AP.Math.Sqr(f);
2415	d2f = -(2fdf);
2416	return;
2417	}
2418	if( k==3 )
2419	{
2420
2421	//
2422	// EX activation function
2423	//
2424	if( (double)(net)>=(double)(0) )
2425	{
2426	net2 = net*net;
2427	arg = net2+1;
2428	root = Math.Sqrt(arg);
2429	f = net+root;
2430	r = net/root;
2431	df = 1+r;
2432	d2f = (root-net*r)/arg;
2433	}
2434	else
2435	{
2436	f = Math.Exp(net);
2437	df = f;
2438	d2f = f;
2439	}
2440	return;
2441	}
2442	if( k==2 )
2443	{
2444	f = Math.Exp(-AP.Math.Sqr(net));
2445	df = -(2netf);
2446	d2f = -(2(f+dfnet));
2447	return;
2448	}
2449	}
2450
2451
2452	/*************************************************************************
2453	Internal subroutine for Hessian calculation.
2454
2455	WARNING!!! Unspeakable math far beyong human capabilities :)
2456	*************************************************************************/
2457	private static void mlphessianbatchinternal(ref multilayerperceptron network,
2458	ref double[,] xy,
2459	int ssize,
2460	bool naturalerr,
2461	ref double e,
2462	ref double[] grad,
2463	ref double[,] h)
2464	{
2465	int nin = 0;
2466	int nout = 0;
2467	int wcount = 0;
2468	int ntotal = 0;
2469	int istart = 0;
2470	int i = 0;
2471	int j = 0;
2472	int k = 0;
2473	int kl = 0;
2474	int offs = 0;
2475	int n1 = 0;
2476	int n2 = 0;
2477	int w1 = 0;
2478	int w2 = 0;
2479	double s = 0;
2480	double t = 0;
2481	double v = 0;
2482	double et = 0;
2483	bool bflag = new bool();
2484	double f = 0;
2485	double df = 0;
2486	double d2f = 0;
2487	double deidyj = 0;
2488	double mx = 0;
2489	double q = 0;
2490	double z = 0;
2491	double s2 = 0;
2492	double expi = 0;
2493	double expj = 0;
2494	double[] x = new double[0];
2495	double[] desiredy = new double[0];
2496	double[] gt = new double[0];
2497	double[] zeros = new double[0];
2498	double[,] rx = new double[0,0];
2499	double[,] ry = new double[0,0];
2500	double[,] rdx = new double[0,0];
2501	double[,] rdy = new double[0,0];
2502	int i_ = 0;
2503	int i1_ = 0;
2504
2505	mlpproperties(ref network, ref nin, ref nout, ref wcount);
2506	ntotal = network.structinfo[3];
2507	istart = network.structinfo[5];
2508
2509	//
2510	// Prepare
2511	//
2512	x = new double[nin-1+1];
2513	desiredy = new double[nout-1+1];
2514	zeros = new double[wcount-1+1];
2515	gt = new double[wcount-1+1];
2516	rx = new double[ntotal+nout-1+1, wcount-1+1];
2517	ry = new double[ntotal+nout-1+1, wcount-1+1];
2518	rdx = new double[ntotal+nout-1+1, wcount-1+1];
2519	rdy = new double[ntotal+nout-1+1, wcount-1+1];
2520	e = 0;
2521	for(i=0; i<=wcount-1; i++)
2522	{
2523	zeros[i] = 0;
2524	}
2525	for(i_=0; i_<=wcount-1;i_++)
2526	{
2527	grad[i_] = zeros[i_];
2528	}
2529	for(i=0; i<=wcount-1; i++)
2530	{
2531	for(i_=0; i_<=wcount-1;i_++)
2532	{
2533	h[i,i_] = zeros[i_];
2534	}
2535	}
2536
2537	//
2538	// Process
2539	//
2540	for(k=0; k<=ssize-1; k++)
2541	{
2542
2543	//
2544	// Process vector with MLPGradN.
2545	// Now Neurons, DFDNET and DError contains results of the last run.
2546	//
2547	for(i_=0; i_<=nin-1;i_++)
2548	{
2549	x[i_] = xy[k,i_];
2550	}
2551	if( mlpissoftmax(ref network) )
2552	{
2553
2554	//
2555	// class labels outputs
2556	//
2557	kl = (int)Math.Round(xy[k,nin]);
2558	for(i=0; i<=nout-1; i++)
2559	{
2560	if( i==kl )
2561	{
2562	desiredy[i] = 1;
2563	}
2564	else
2565	{
2566	desiredy[i] = 0;
2567	}
2568	}
2569	}
2570	else
2571	{
2572
2573	//
2574	// real outputs
2575	//
2576	i1_ = (nin) - (0);
2577	for(i_=0; i_<=nout-1;i_++)
2578	{
2579	desiredy[i_] = xy[k,i_+i1_];
2580	}
2581	}
2582	if( naturalerr )
2583	{
2584	mlpgradn(ref network, ref x, ref desiredy, ref et, ref gt);
2585	}
2586	else
2587	{
2588	mlpgrad(ref network, ref x, ref desiredy, ref et, ref gt);
2589	}
2590
2591	//
2592	// grad, error
2593	//
2594	e = e+et;
2595	for(i_=0; i_<=wcount-1;i_++)
2596	{
2597	grad[i_] = grad[i_] + gt[i_];
2598	}
2599
2600	//
2601	// Hessian.
2602	// Forward pass of the R-algorithm
2603	//
2604	for(i=0; i<=ntotal-1; i++)
2605	{
2606	offs = istart+i*nfieldwidth;
2607	for(i_=0; i_<=wcount-1;i_++)
2608	{
2609	rx[i,i_] = zeros[i_];
2610	}
2611	for(i_=0; i_<=wcount-1;i_++)
2612	{
2613	ry[i,i_] = zeros[i_];
2614	}
2615	if( network.structinfo[offs+0]>0 )
2616	{
2617
2618	//
2619	// Activation function
2620	//
2621	n1 = network.structinfo[offs+2];
2622	for(i_=0; i_<=wcount-1;i_++)
2623	{
2624	rx[i,i_] = ry[n1,i_];
2625	}
2626	v = network.dfdnet[i];
2627	for(i_=0; i_<=wcount-1;i_++)
2628	{
2629	ry[i,i_] = v*rx[i,i_];
2630	}
2631	}
2632	if( network.structinfo[offs+0]==0 )
2633	{
2634
2635	//
2636	// Adaptive summator
2637	//
2638	n1 = network.structinfo[offs+2];
2639	n2 = n1+network.structinfo[offs+1]-1;
2640	w1 = network.structinfo[offs+3];
2641	w2 = w1+network.structinfo[offs+1]-1;
2642	for(j=n1; j<=n2; j++)
2643	{
2644	v = network.weights[w1+j-n1];
2645	for(i_=0; i_<=wcount-1;i_++)
2646	{
2647	rx[i,i_] = rx[i,i_] + v*ry[j,i_];
2648	}
2649	rx[i,w1+j-n1] = rx[i,w1+j-n1]+network.neurons[j];
2650	}
2651	for(i_=0; i_<=wcount-1;i_++)
2652	{
2653	ry[i,i_] = rx[i,i_];
2654	}
2655	}
2656	if( network.structinfo[offs+0]<0 )
2657	{
2658	bflag = true;
2659	if( network.structinfo[offs+0]==-2 )
2660	{
2661
2662	//
2663	// input neuron, left unchanged
2664	//
2665	bflag = false;
2666	}
2667	if( network.structinfo[offs+0]==-3 )
2668	{
2669
2670	//
2671	// "-1" neuron, left unchanged
2672	//
2673	bflag = false;
2674	}
2675	if( network.structinfo[offs+0]==-4 )
2676	{
2677
2678	//
2679	// "0" neuron, left unchanged
2680	//
2681	bflag = false;
2682	}
2683	System.Diagnostics.Debug.Assert(!bflag, "MLPHessianNBatch: internal error - unknown neuron type!");
2684	}
2685	}
2686
2687	//
2688	// Hessian. Backward pass of the R-algorithm.
2689	//
2690	// Stage 1. Initialize RDY
2691	//
2692	for(i=0; i<=ntotal+nout-1; i++)
2693	{
2694	for(i_=0; i_<=wcount-1;i_++)
2695	{
2696	rdy[i,i_] = zeros[i_];
2697	}
2698	}
2699	if( network.structinfo[6]==0 )
2700	{
2701
2702	//
2703	// Standardisation.
2704	//
2705	// In context of the Hessian calculation standardisation
2706	// is considered as additional layer with weightless
2707	// activation function:
2708	//
2709	// F(NET) := Sigma*NET
2710	//
2711	// So we add one more layer to forward pass, and
2712	// make forward/backward pass through this layer.
2713	//
2714	for(i=0; i<=nout-1; i++)
2715	{
2716	n1 = ntotal-nout+i;
2717	n2 = ntotal+i;
2718
2719	//
2720	// Forward pass from N1 to N2
2721	//
2722	for(i_=0; i_<=wcount-1;i_++)
2723	{
2724	rx[n2,i_] = ry[n1,i_];
2725	}
2726	v = network.columnsigmas[nin+i];
2727	for(i_=0; i_<=wcount-1;i_++)
2728	{
2729	ry[n2,i_] = v*rx[n2,i_];
2730	}
2731
2732	//
2733	// Initialization of RDY
2734	//
2735	for(i_=0; i_<=wcount-1;i_++)
2736	{
2737	rdy[n2,i_] = ry[n2,i_];
2738	}
2739
2740	//
2741	// Backward pass from N2 to N1:
2742	// 1. Calculate R(dE/dX).
2743	// 2. No R(dE/dWij) is needed since weight of activation neuron
2744	// is fixed to 1. So we can update R(dE/dY) for
2745	// the connected neuron (note that Vij=0, Wij=1)
2746	//
2747	df = network.columnsigmas[nin+i];
2748	for(i_=0; i_<=wcount-1;i_++)
2749	{
2750	rdx[n2,i_] = df*rdy[n2,i_];
2751	}
2752	for(i_=0; i_<=wcount-1;i_++)
2753	{
2754	rdy[n1,i_] = rdy[n1,i_] + rdx[n2,i_];
2755	}
2756	}
2757	}
2758	else
2759	{
2760
2761	//
2762	// Softmax.
2763	//
2764	// Initialize RDY using generalized expression for ei'(yi)
2765	// (see expression (9) from p. 5 of "Fast Exact Multiplication by the Hessian").
2766	//
2767	// When we are working with softmax network, generalized
2768	// expression for ei'(yi) is used because softmax
2769	// normalization leads to ei, which depends on all y's
2770	//
2771	if( naturalerr )
2772	{
2773
2774	//
2775	// softmax + cross-entropy.
2776	// We have:
2777	//
2778	// S = sum(exp(yk)),
2779	// ei = sum(trn)*exp(yi)/S-trn_i
2780	//
2781	// j=i: d(ei)/d(yj) = Texp(yi)(S-exp(yi))/S^2
2782	// j<>i: d(ei)/d(yj) = -Texp(yi)exp(yj)/S^2
2783	//
2784	t = 0;
2785	for(i=0; i<=nout-1; i++)
2786	{
2787	t = t+desiredy[i];
2788	}
2789	mx = network.neurons[ntotal-nout];
2790	for(i=0; i<=nout-1; i++)
2791	{
2792	mx = Math.Max(mx, network.neurons[ntotal-nout+i]);
2793	}
2794	s = 0;
2795	for(i=0; i<=nout-1; i++)
2796	{
2797	network.nwbuf[i] = Math.Exp(network.neurons[ntotal-nout+i]-mx);
2798	s = s+network.nwbuf[i];
2799	}
2800	for(i=0; i<=nout-1; i++)
2801	{
2802	for(j=0; j<=nout-1; j++)
2803	{
2804	if( j==i )
2805	{
2806	deidyj = tnetwork.nwbuf[i](s-network.nwbuf[i])/AP.Math.Sqr(s);
2807	for(i_=0; i_<=wcount-1;i_++)
2808	{
2809	rdy[ntotal-nout+i,i_] = rdy[ntotal-nout+i,i_] + deidyj*ry[ntotal-nout+i,i_];
2810	}
2811	}
2812	else
2813	{
2814	deidyj = -(tnetwork.nwbuf[i]network.nwbuf[j]/AP.Math.Sqr(s));
2815	for(i_=0; i_<=wcount-1;i_++)
2816	{
2817	rdy[ntotal-nout+i,i_] = rdy[ntotal-nout+i,i_] + deidyj*ry[ntotal-nout+j,i_];
2818	}
2819	}
2820	}
2821	}
2822	}
2823	else
2824	{
2825
2826	//
2827	// For a softmax + squared error we have expression
2828	// far beyond human imagination so we dont even try
2829	// to comment on it. Just enjoy the code...
2830	//
2831	// P.S. That's why "natural error" is called "natural" -
2832	// compact beatiful expressions, fast code....
2833	//
2834	mx = network.neurons[ntotal-nout];
2835	for(i=0; i<=nout-1; i++)
2836	{
2837	mx = Math.Max(mx, network.neurons[ntotal-nout+i]);
2838	}
2839	s = 0;
2840	s2 = 0;
2841	for(i=0; i<=nout-1; i++)
2842	{
2843	network.nwbuf[i] = Math.Exp(network.neurons[ntotal-nout+i]-mx);
2844	s = s+network.nwbuf[i];
2845	s2 = s2+AP.Math.Sqr(network.nwbuf[i]);
2846	}
2847	q = 0;
2848	for(i=0; i<=nout-1; i++)
2849	{
2850	q = q+(network.y[i]-desiredy[i])*network.nwbuf[i];
2851	}
2852	for(i=0; i<=nout-1; i++)
2853	{
2854	z = -q+(network.y[i]-desiredy[i])*s;
2855	expi = network.nwbuf[i];
2856	for(j=0; j<=nout-1; j++)
2857	{
2858	expj = network.nwbuf[j];
2859	if( j==i )
2860	{
2861	deidyj = expi/AP.Math.Sqr(s)((z+expi)(s-2expi)/s+expis2/AP.Math.Sqr(s));
2862	}
2863	else
2864	{
2865	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]));
2866	}
2867	for(i_=0; i_<=wcount-1;i_++)
2868	{
2869	rdy[ntotal-nout+i,i_] = rdy[ntotal-nout+i,i_] + deidyj*ry[ntotal-nout+j,i_];
2870	}
2871	}
2872	}
2873	}
2874	}
2875
2876	//
2877	// Hessian. Backward pass of the R-algorithm
2878	//
2879	// Stage 2. Process.
2880	//
2881	for(i=ntotal-1; i>=0; i--)
2882	{
2883
2884	//
2885	// Possible variants:
2886	// 1. Activation function
2887	// 2. Adaptive summator
2888	// 3. Special neuron
2889	//
2890	offs = istart+i*nfieldwidth;
2891	if( network.structinfo[offs+0]>0 )
2892	{
2893	n1 = network.structinfo[offs+2];
2894
2895	//
2896	// First, calculate R(dE/dX).
2897	//
2898	mlpactivationfunction(network.neurons[n1], network.structinfo[offs+0], ref f, ref df, ref d2f);
2899	v = d2f*network.derror[i];
2900	for(i_=0; i_<=wcount-1;i_++)
2901	{
2902	rdx[i,i_] = df*rdy[i,i_];
2903	}
2904	for(i_=0; i_<=wcount-1;i_++)
2905	{
2906	rdx[i,i_] = rdx[i,i_] + v*rx[i,i_];
2907	}
2908
2909	//
2910	// No R(dE/dWij) is needed since weight of activation neuron
2911	// is fixed to 1.
2912	//
2913	// So we can update R(dE/dY) for the connected neuron.
2914	// (note that Vij=0, Wij=1)
2915	//
2916	for(i_=0; i_<=wcount-1;i_++)
2917	{
2918	rdy[n1,i_] = rdy[n1,i_] + rdx[i,i_];
2919	}
2920	}
2921	if( network.structinfo[offs+0]==0 )
2922	{
2923
2924	//
2925	// Adaptive summator
2926	//
2927	n1 = network.structinfo[offs+2];
2928	n2 = n1+network.structinfo[offs+1]-1;
2929	w1 = network.structinfo[offs+3];
2930	w2 = w1+network.structinfo[offs+1]-1;
2931
2932	//
2933	// First, calculate R(dE/dX).
2934	//
2935	for(i_=0; i_<=wcount-1;i_++)
2936	{
2937	rdx[i,i_] = rdy[i,i_];
2938	}
2939
2940	//
2941	// Then, calculate R(dE/dWij)
2942	//
2943	for(j=w1; j<=w2; j++)
2944	{
2945	v = network.neurons[n1+j-w1];
2946	for(i_=0; i_<=wcount-1;i_++)
2947	{
2948	h[j,i_] = h[j,i_] + v*rdx[i,i_];
2949	}
2950	v = network.derror[i];
2951	for(i_=0; i_<=wcount-1;i_++)
2952	{
2953	h[j,i_] = h[j,i_] + v*ry[n1+j-w1,i_];
2954	}
2955	}
2956
2957	//
2958	// And finally, update R(dE/dY) for connected neurons.
2959	//
2960	for(j=w1; j<=w2; j++)
2961	{
2962	v = network.weights[j];
2963	for(i_=0; i_<=wcount-1;i_++)
2964	{
2965	rdy[n1+j-w1,i_] = rdy[n1+j-w1,i_] + v*rdx[i,i_];
2966	}
2967	rdy[n1+j-w1,j] = rdy[n1+j-w1,j]+network.derror[i];
2968	}
2969	}
2970	if( network.structinfo[offs+0]<0 )
2971	{
2972	bflag = false;
2973	if( network.structinfo[offs+0]==-2 \| network.structinfo[offs+0]==-3 \| network.structinfo[offs+0]==-4 )
2974	{
2975
2976	//
2977	// Special neuron type, no back-propagation required
2978	//
2979	bflag = true;
2980	}
2981	System.Diagnostics.Debug.Assert(bflag, "MLPHessianNBatch: unknown neuron type!");
2982	}
2983	}
2984	}
2985	}
2986
2987
2988	/*************************************************************************
2989	Internal subroutine
2990
2991	Network must be processed by MLPProcess on X
2992	*************************************************************************/
2993	private static void mlpinternalcalculategradient(ref multilayerperceptron network,
2994	ref double[] neurons,
2995	ref double[] weights,
2996	ref double[] derror,
2997	ref double[] grad,
2998	bool naturalerrorfunc)
2999	{
3000	int i = 0;
3001	int n1 = 0;
3002	int n2 = 0;
3003	int w1 = 0;
3004	int w2 = 0;
3005	int ntotal = 0;
3006	int istart = 0;
3007	int nin = 0;
3008	int nout = 0;
3009	int offs = 0;
3010	double dedf = 0;
3011	double dfdnet = 0;
3012	double v = 0;
3013	double fown = 0;
3014	double deown = 0;
3015	double net = 0;
3016	double mx = 0;
3017	bool bflag = new bool();
3018	int i_ = 0;
3019	int i1_ = 0;
3020
3021
3022	//
3023	// Read network geometry
3024	//
3025	nin = network.structinfo[1];
3026	nout = network.structinfo[2];
3027	ntotal = network.structinfo[3];
3028	istart = network.structinfo[5];
3029
3030	//
3031	// Pre-processing of dError/dOut:
3032	// from dError/dOut(normalized) to dError/dOut(non-normalized)
3033	//
3034	System.Diagnostics.Debug.Assert(network.structinfo[6]==0 \| network.structinfo[6]==1, "MLPInternalCalculateGradient: unknown normalization type!");
3035	if( network.structinfo[6]==1 )
3036	{
3037
3038	//
3039	// Softmax
3040	//
3041	if( !naturalerrorfunc )
3042	{
3043	mx = network.neurons[ntotal-nout];
3044	for(i=0; i<=nout-1; i++)
3045	{
3046	mx = Math.Max(mx, network.neurons[ntotal-nout+i]);
3047	}
3048	net = 0;
3049	for(i=0; i<=nout-1; i++)
3050	{
3051	network.nwbuf[i] = Math.Exp(network.neurons[ntotal-nout+i]-mx);
3052	net = net+network.nwbuf[i];
3053	}
3054	i1_ = (0)-(ntotal-nout);
3055	v = 0.0;
3056	for(i_=ntotal-nout; i_<=ntotal-1;i_++)
3057	{
3058	v += network.derror[i_]*network.nwbuf[i_+i1_];
3059	}
3060	for(i=0; i<=nout-1; i++)
3061	{
3062	fown = network.nwbuf[i];
3063	deown = network.derror[ntotal-nout+i];
3064	network.nwbuf[nout+i] = (-v+deownfown+deown(net-fown))*fown/AP.Math.Sqr(net);
3065	}
3066	for(i=0; i<=nout-1; i++)
3067	{
3068	network.derror[ntotal-nout+i] = network.nwbuf[nout+i];
3069	}
3070	}
3071	}
3072	else
3073	{
3074
3075	//
3076	// Un-standardisation
3077	//
3078	for(i=0; i<=nout-1; i++)
3079	{
3080	network.derror[ntotal-nout+i] = network.derror[ntotal-nout+i]*network.columnsigmas[nin+i];
3081	}
3082	}
3083
3084	//
3085	// Backpropagation
3086	//
3087	for(i=ntotal-1; i>=0; i--)
3088	{
3089
3090	//
3091	// Extract info
3092	//
3093	offs = istart+i*nfieldwidth;
3094	if( network.structinfo[offs+0]>0 )
3095	{
3096
3097	//
3098	// Activation function
3099	//
3100	dedf = network.derror[i];
3101	dfdnet = network.dfdnet[i];
3102	derror[network.structinfo[offs+2]] = derror[network.structinfo[offs+2]]+dedf*dfdnet;
3103	}
3104	if( network.structinfo[offs+0]==0 )
3105	{
3106
3107	//
3108	// Adaptive summator
3109	//
3110	n1 = network.structinfo[offs+2];
3111	n2 = n1+network.structinfo[offs+1]-1;
3112	w1 = network.structinfo[offs+3];
3113	w2 = w1+network.structinfo[offs+1]-1;
3114	dedf = network.derror[i];
3115	dfdnet = 1.0;
3116	v = dedf*dfdnet;
3117	i1_ = (n1) - (w1);
3118	for(i_=w1; i_<=w2;i_++)
3119	{
3120	grad[i_] = v*neurons[i_+i1_];
3121	}
3122	i1_ = (w1) - (n1);
3123	for(i_=n1; i_<=n2;i_++)
3124	{
3125	derror[i_] = derror[i_] + v*weights[i_+i1_];
3126	}
3127	}
3128	if( network.structinfo[offs+0]<0 )
3129	{
3130	bflag = false;
3131	if( network.structinfo[offs+0]==-2 \| network.structinfo[offs+0]==-3 \| network.structinfo[offs+0]==-4 )
3132	{
3133
3134	//
3135	// Special neuron type, no back-propagation required
3136	//
3137	bflag = true;
3138	}
3139	System.Diagnostics.Debug.Assert(bflag, "MLPInternalCalculateGradient: unknown neuron type!");
3140	}
3141	}
3142	}
3143
3144
3145	/*************************************************************************
3146	Internal subroutine, chunked gradient
3147	*************************************************************************/
3148	private static void mlpchunkedgradient(ref multilayerperceptron network,
3149	ref double[,] xy,
3150	int cstart,
3151	int csize,
3152	ref double e,
3153	ref double[] grad,
3154	bool naturalerrorfunc)
3155	{
3156	int i = 0;
3157	int j = 0;
3158	int k = 0;
3159	int kl = 0;
3160	int n1 = 0;
3161	int n2 = 0;
3162	int w1 = 0;
3163	int w2 = 0;
3164	int c1 = 0;
3165	int c2 = 0;
3166	int ntotal = 0;
3167	int nin = 0;
3168	int nout = 0;
3169	int offs = 0;
3170	double f = 0;
3171	double df = 0;
3172	double d2f = 0;
3173	double v = 0;
3174	double s = 0;
3175	double fown = 0;
3176	double deown = 0;
3177	double net = 0;
3178	double lnnet = 0;
3179	double mx = 0;
3180	bool bflag = new bool();
3181	int istart = 0;
3182	int ineurons = 0;
3183	int idfdnet = 0;
3184	int iderror = 0;
3185	int izeros = 0;
3186	int i_ = 0;
3187	int i1_ = 0;
3188
3189
3190	//
3191	// Read network geometry, prepare data
3192	//
3193	nin = network.structinfo[1];
3194	nout = network.structinfo[2];
3195	ntotal = network.structinfo[3];
3196	istart = network.structinfo[5];
3197	c1 = cstart;
3198	c2 = cstart+csize-1;
3199	ineurons = 0;
3200	idfdnet = ntotal;
3201	iderror = 2*ntotal;
3202	izeros = 3*ntotal;
3203	for(j=0; j<=csize-1; j++)
3204	{
3205	network.chunks[izeros,j] = 0;
3206	}
3207
3208	//
3209	// Forward pass:
3210	// 1. Load inputs from XY to Chunks[0:NIn-1,0:CSize-1]
3211	// 2. Forward pass
3212	//
3213	for(i=0; i<=nin-1; i++)
3214	{
3215	for(j=0; j<=csize-1; j++)
3216	{
3217	if( (double)(network.columnsigmas[i])!=(double)(0) )
3218	{
3219	network.chunks[i,j] = (xy[c1+j,i]-network.columnmeans[i])/network.columnsigmas[i];
3220	}
3221	else
3222	{
3223	network.chunks[i,j] = xy[c1+j,i]-network.columnmeans[i];
3224	}
3225	}
3226	}
3227	for(i=0; i<=ntotal-1; i++)
3228	{
3229	offs = istart+i*nfieldwidth;
3230	if( network.structinfo[offs+0]>0 )
3231	{
3232
3233	//
3234	// Activation function:
3235	// * calculate F vector, F(i) = F(NET(i))
3236	//
3237	n1 = network.structinfo[offs+2];
3238	for(i_=0; i_<=csize-1;i_++)
3239	{
3240	network.chunks[i,i_] = network.chunks[n1,i_];
3241	}
3242	for(j=0; j<=csize-1; j++)
3243	{
3244	mlpactivationfunction(network.chunks[i,j], network.structinfo[offs+0], ref f, ref df, ref d2f);
3245	network.chunks[i,j] = f;
3246	network.chunks[idfdnet+i,j] = df;
3247	}
3248	}
3249	if( network.structinfo[offs+0]==0 )
3250	{
3251
3252	//
3253	// Adaptive summator:
3254	// * calculate NET vector, NET(i) = SUM(W(j,i)*Neurons(j),j=N1..N2)
3255	//
3256	n1 = network.structinfo[offs+2];
3257	n2 = n1+network.structinfo[offs+1]-1;
3258	w1 = network.structinfo[offs+3];
3259	w2 = w1+network.structinfo[offs+1]-1;
3260	for(i_=0; i_<=csize-1;i_++)
3261	{
3262	network.chunks[i,i_] = network.chunks[izeros,i_];
3263	}
3264	for(j=n1; j<=n2; j++)
3265	{
3266	v = network.weights[w1+j-n1];
3267	for(i_=0; i_<=csize-1;i_++)
3268	{
3269	network.chunks[i,i_] = network.chunks[i,i_] + v*network.chunks[j,i_];
3270	}
3271	}
3272	}
3273	if( network.structinfo[offs+0]<0 )
3274	{
3275	bflag = false;
3276	if( network.structinfo[offs+0]==-2 )
3277	{
3278
3279	//
3280	// input neuron, left unchanged
3281	//
3282	bflag = true;
3283	}
3284	if( network.structinfo[offs+0]==-3 )
3285	{
3286
3287	//
3288	// "-1" neuron
3289	//
3290	for(k=0; k<=csize-1; k++)
3291	{
3292	network.chunks[i,k] = -1;
3293	}
3294	bflag = true;
3295	}
3296	if( network.structinfo[offs+0]==-4 )
3297	{
3298
3299	//
3300	// "0" neuron
3301	//
3302	for(k=0; k<=csize-1; k++)
3303	{
3304	network.chunks[i,k] = 0;
3305	}
3306	bflag = true;
3307	}
3308	System.Diagnostics.Debug.Assert(bflag, "MLPChunkedGradient: internal error - unknown neuron type!");
3309	}
3310	}
3311
3312	//
3313	// Post-processing, error, dError/dOut
3314	//
3315	for(i=0; i<=ntotal-1; i++)
3316	{
3317	for(i_=0; i_<=csize-1;i_++)
3318	{
3319	network.chunks[iderror+i,i_] = network.chunks[izeros,i_];
3320	}
3321	}
3322	System.Diagnostics.Debug.Assert(network.structinfo[6]==0 \| network.structinfo[6]==1, "MLPChunkedGradient: unknown normalization type!");
3323	if( network.structinfo[6]==1 )
3324	{
3325
3326	//
3327	// Softmax output, classification network.
3328	//
3329	// For each K = 0..CSize-1 do:
3330	// 1. place exp(outputs[k]) to NWBuf[0:NOut-1]
3331	// 2. place sum(exp(..)) to NET
3332	// 3. calculate dError/dOut and place it to the second block of Chunks
3333	//
3334	for(k=0; k<=csize-1; k++)
3335	{
3336
3337	//
3338	// Normalize
3339	//
3340	mx = network.chunks[ntotal-nout,k];
3341	for(i=1; i<=nout-1; i++)
3342	{
3343	mx = Math.Max(mx, network.chunks[ntotal-nout+i,k]);
3344	}
3345	net = 0;
3346	for(i=0; i<=nout-1; i++)
3347	{
3348	network.nwbuf[i] = Math.Exp(network.chunks[ntotal-nout+i,k]-mx);
3349	net = net+network.nwbuf[i];
3350	}
3351
3352	//
3353	// Calculate error function and dError/dOut
3354	//
3355	if( naturalerrorfunc )
3356	{
3357
3358	//
3359	// Natural error func.
3360	//
3361	//
3362	s = 1;
3363	lnnet = Math.Log(net);
3364	kl = (int)Math.Round(xy[cstart+k,nin]);
3365	for(i=0; i<=nout-1; i++)
3366	{
3367	if( i==kl )
3368	{
3369	v = 1;
3370	}
3371	else
3372	{
3373	v = 0;
3374	}
3375	network.chunks[iderror+ntotal-nout+i,k] = s*network.nwbuf[i]/net-v;
3376	e = e+safecrossentropy(v, network.nwbuf[i]/net);
3377	}
3378	}
3379	else
3380	{
3381
3382	//
3383	// Least squares error func
3384	// Error, dError/dOut(normalized)
3385	//
3386	kl = (int)Math.Round(xy[cstart+k,nin]);
3387	for(i=0; i<=nout-1; i++)
3388	{
3389	if( i==kl )
3390	{
3391	v = network.nwbuf[i]/net-1;
3392	}
3393	else
3394	{
3395	v = network.nwbuf[i]/net;
3396	}
3397	network.nwbuf[nout+i] = v;
3398	e = e+AP.Math.Sqr(v)/2;
3399	}
3400
3401	//
3402	// From dError/dOut(normalized) to dError/dOut(non-normalized)
3403	//
3404	i1_ = (0)-(nout);
3405	v = 0.0;
3406	for(i_=nout; i_<=2*nout-1;i_++)
3407	{
3408	v += network.nwbuf[i_]*network.nwbuf[i_+i1_];
3409	}
3410	for(i=0; i<=nout-1; i++)
3411	{
3412	fown = network.nwbuf[i];
3413	deown = network.nwbuf[nout+i];
3414	network.chunks[iderror+ntotal-nout+i,k] = (-v+deownfown+deown(net-fown))*fown/AP.Math.Sqr(net);
3415	}
3416	}
3417	}
3418	}
3419	else
3420	{
3421
3422	//
3423	// Normal output, regression network
3424	//
3425	// For each K = 0..CSize-1 do:
3426	// 1. calculate dError/dOut and place it to the second block of Chunks
3427	//
3428	for(i=0; i<=nout-1; i++)
3429	{
3430	for(j=0; j<=csize-1; j++)
3431	{
3432	v = network.chunks[ntotal-nout+i,j]*network.columnsigmas[nin+i]+network.columnmeans[nin+i]-xy[cstart+j,nin+i];
3433	network.chunks[iderror+ntotal-nout+i,j] = v*network.columnsigmas[nin+i];
3434	e = e+AP.Math.Sqr(v)/2;
3435	}
3436	}
3437	}
3438
3439	//
3440	// Backpropagation
3441	//
3442	for(i=ntotal-1; i>=0; i--)
3443	{
3444
3445	//
3446	// Extract info
3447	//
3448	offs = istart+i*nfieldwidth;
3449	if( network.structinfo[offs+0]>0 )
3450	{
3451
3452	//
3453	// Activation function
3454	//
3455	n1 = network.structinfo[offs+2];
3456	for(k=0; k<=csize-1; k++)
3457	{
3458	network.chunks[iderror+i,k] = network.chunks[iderror+i,k]*network.chunks[idfdnet+i,k];
3459	}
3460	for(i_=0; i_<=csize-1;i_++)
3461	{
3462	network.chunks[iderror+n1,i_] = network.chunks[iderror+n1,i_] + network.chunks[iderror+i,i_];
3463	}
3464	}
3465	if( network.structinfo[offs+0]==0 )
3466	{
3467
3468	//
3469	// "Normal" activation function
3470	//
3471	n1 = network.structinfo[offs+2];
3472	n2 = n1+network.structinfo[offs+1]-1;
3473	w1 = network.structinfo[offs+3];
3474	w2 = w1+network.structinfo[offs+1]-1;
3475	for(j=w1; j<=w2; j++)
3476	{
3477	v = 0.0;
3478	for(i_=0; i_<=csize-1;i_++)
3479	{
3480	v += network.chunks[n1+j-w1,i_]*network.chunks[iderror+i,i_];
3481	}
3482	grad[j] = grad[j]+v;
3483	}
3484	for(j=n1; j<=n2; j++)
3485	{
3486	v = network.weights[w1+j-n1];
3487	for(i_=0; i_<=csize-1;i_++)
3488	{
3489	network.chunks[iderror+j,i_] = network.chunks[iderror+j,i_] + v*network.chunks[iderror+i,i_];
3490	}
3491	}
3492	}
3493	if( network.structinfo[offs+0]<0 )
3494	{
3495	bflag = false;
3496	if( network.structinfo[offs+0]==-2 \| network.structinfo[offs+0]==-3 \| network.structinfo[offs+0]==-4 )
3497	{
3498
3499	//
3500	// Special neuron type, no back-propagation required
3501	//
3502	bflag = true;
3503	}
3504	System.Diagnostics.Debug.Assert(bflag, "MLPInternalCalculateGradient: unknown neuron type!");
3505	}
3506	}
3507	}
3508
3509
3510	/*************************************************************************
3511	Returns T*Ln(T/Z), guarded against overflow/underflow.
3512	Internal subroutine.
3513	*************************************************************************/
3514	private static double safecrossentropy(double t,
3515	double z)
3516	{
3517	double result = 0;
3518	double r = 0;
3519
3520	if( (double)(t)==(double)(0) )
3521	{
3522	result = 0;
3523	}
3524	else
3525	{
3526	if( (double)(Math.Abs(z))>(double)(1) )
3527	{
3528
3529	//
3530	// Shouldn't be the case with softmax,
3531	// but we just want to be sure.
3532	//
3533	if( (double)(t/z)==(double)(0) )
3534	{
3535	r = AP.Math.MinRealNumber;
3536	}
3537	else
3538	{
3539	r = t/z;
3540	}
3541	}
3542	else
3543	{
3544
3545	//
3546	// Normal case
3547	//
3548	if( (double)(z)==(double)(0) \| (double)(Math.Abs(t))>=(double)(AP.Math.MaxRealNumber*Math.Abs(z)) )
3549	{
3550	r = AP.Math.MaxRealNumber;
3551	}
3552	else
3553	{
3554	r = t/z;
3555	}
3556	}
3557	result = t*Math.Log(r);
3558	}
3559	return result;
3560	}
3561	}
3562	}

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

Download in other formats:

Update cookies preferences