
# This is a shell archive.  Remove anything before this line,
# then unpack it by saving it in a file and typing "sh file".
#
# Wrapped by Dr. Philip Goodman <goodman@unr.edu> on Wed May 26 20:57:54 1993
#
# This archive contains:
#	NevProp/NP_READ.ME	NevProp/Makefile	
#	NevProp/np.c		NevProp/cind.c		
#	NevProp/create.c	NevProp/rank.c		
#	NevProp/cind.h		NevProp/const.h		
#	NevProp/create.h	NevProp/do.h		
#	NevProp/rank.h		NevProp/protoTyp.h	
#	NevProp/iris.net	
#

LANG=""; export LANG
PATH=/bin:/usr/bin:$PATH; export PATH

echo x - Made NevProp directory
mkdir NevProp

echo x - NevProp/NP_READ.ME
sed 's/^@//' >NevProp/NP_READ.ME <<'@EOF'
 /* * * * * * * This READ.ME describes NevProp, version 1.xx * * * * * * 
                    - Last updated May 23, 1993 -
  
 .- - - - QUICKSTART - - FOR THOSE WHO DON'T LIKE TO READ MANUALS - - - .
 |                                                                      |
 |  UNIX:                                                               |
 |    1. ftp the 'npxxx.shar' file  (xxx depends on current version)    |
 |    2. type 'sh npxxx.shar'  (creates 'NevProp' directory - cd to it) |
 |    3. type 'make'  (should create an executable 'np' program)        |
 |    4. type 'np'                                                      |
 |    5. at prompt for .net filename, type 'iris'  (supplied)           |
 |    6. explore!                                                       |
 |                                                                      |
 |  MAC & PCs:                                                          |
 |    1. execute the program  (assumes you got it somehow)              |
 |    2. at prompt for .net filename type 'iris'  (supplied)            |
 |    3. explore!                                                       |
 ._ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ .

   NevProp is a general purpose back-propagation program written in C for
   UNIX, Macintosh, and DOS. The original version was Quickprop 1.0 by
   Scott Fahlman, as translated from Common Lisp into C by Terry Regier.

   The quickprop algorithm itself was not substantively changed, but we 
   inserted options to force gradient descent (per-epoch or per-pattern)
   and added generalization & stopped training, c index, and interface
   enhancements.
  
   Questions and comments regarding the Quickprop algorithm itself should 
   be directed to Scott Fahlman at sef+@cs.cmu.edu. Issues pertaining 
   to NevProp code or functionality should come our way. 
  
    Phil Goodman (goodman@unr.edu)
    David Rosen   
    Allen Plummer
      University of Nevada Center for Biomedical Modeling Research
      Washoe Medical Center, H-166, 77 Pringle Way, Reno, NV  89520 
      702-328-4867  FAX: 328-4111
  
    We express our appreciation for their support to Medicine Chairman
    Stan Shane, Medical School Dean Bob Daugherty, Graduate School Dean
    Ken Hunter, and Washoe Medical Center CEO Bob Burn.
  
    And thanks to Dwight Egbert for the DOS compilation.

************************************************************************

  CONTENTS OF "NP_READ.ME":
 
   &&A. AVAILABILITY & BUG REPORTS
   &&B. RATIONALE for NevProp
   &&C. FEATURES of NevProp
   &&D. WHAT IS THE c index?
   &&E. COMPILING NevProp / Makefile
   &&F. .net FILE and PARAMETER DESCRIPTIONS
   &&G. SAMPLE OUTPUT OF RUNS on "iris.net"
   &&H. SAMPLE CONTENTS of .wts and .prd saved files
  
************************************************************************
 
  &&A. AVAILABILITY & BUG REPORTS...
 
   . - - - - - - - - - - - - BUG REPORTS - - - - - - - - - - - - - - .
   |                                                                 |
   |  Before reporting bugs, please check the unr machine for the    |
   |  "nevprop.bug" file for bug reports, fixes, and upgrade plans.  |
   |                                                                 |
   . _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ .


  3 SOURCES FOR NevProp:

  1. The most updated version of NevProp will be made available by
     anonymous ftp from the University of Nevada, Reno:
 
    a. Create an FTP connection from wherever you are to machine
       "unssun.scs.unr.edu".  The internet address of this machine
       is 134.197.10.128 for those who need it.
 
    b. Log in as user "anonymous" with your own ID as password.  
 
    c. Change remote directory to "pub/goodman/nevpropdir".
 
    d. At this point FTP should be able to get a listing of files
       (The exact FTP commands depend on your local FTP server.)
 
 
  2. NevProp and the original Quickprop are available as C code
     by anonymous ftp from CMU: 
 
    a. Create an FTP connection from wherever you are to machine
       "ftp.cs.cmu.edu".  The internet address of this machine is
       128.2.206.173, for those who need it.
 
    b. Log in as user "anonymous" with your own ID as password.  
       You may see an error message that says "filenames may not
       have /.. in them" or something like that.  Just ignore it.
 
    c. Change remote directory to "/afs/cs/project/connect/code".
       NOTE: You must do this in a single operation.  Some of the
       super directories on this path are protected against outside
       users.
 
    d. At this point FTP should be able to get a listing of files
       in this directory with "dir" & fetch the ones you want with "get".
       (The exact FTP commands depend on your local FTP server.)
 
 
  3. MACINTOSH and DOS versions of NevProp 1.xx:
  
    May be available in the unr anonymous ftp above. If you are unable
    to transfer the files, just send a formatted 3.5 inch floppy disk
    (specifying Mac or DOS) with a self-addressed stamped return
    envelope to Phil Goodman at the CBMR address above, including
    a brief description of your research.
 
************************************************************************
 
  &&B. RATIONALE for NevProp
 
  NevProp was developed at the Univ. of Nevada Center for Biomedical 
  Modeling Research for our campus-wide Neural Network Interest Group.
  Because of the diverse interests and expertise of the Interest Group
  members, NevProp was designed to be user-friendly and operate on 
  multiple platforms, yet be sufficiently powerful to run massive data 
  sets and automatically monitor the training process for over-fitting 
  (i.e., to improve external generalization).
  
  We started with the code carrying QuickProp v1.0 by Scott Fahlman,
  as translated into C by Terry Regier. We then added the interface
  and functional features described below.
  
************************************************************************
 
  &&C. FEATURES of NevProp version 1
   
o UNLIMITED (except by machine memory) number of input PATTERNS;

o UNLIMITED number of input, hidden, and output UNITS;

o Arbitrary CONNECTIONS among the various layers' units;

o Clock-time or user-specified RANDOM SEED for initial random weights;

o Choice of regular GRADIENT DESCENT or QUICKPROP;

o Choice of LOGISTIC or TANH activation functions;

o Choice of PER-EPOCH or PER-PATTERN (stochastic) weight updating;

o GENERALIZATION to a test dataset;

o AUTOMATICALLY STOPPED TRAINING based on generalization;

o RETENTION of best-generalizing weights and predictions;

o Simple but useful bar GRAPH to show smoothness of generalization;

o SAVING of results to a file while working interactively;

o SAVING of weights file and reloading for continued training;

o PREDICTION-only on datasets by applying an existing weights file;

o In addition to RMS error, the concordance, or c index is displayed.
  The c index shows the correctness of the RELATIVE ordering
  of predictions AMONG the cases; ie, it considers all possible PAIRS
  of vectors. This statistic is identical to the area under the
  receiver operating characteristic (ROC) curve, widely used in
  technology assessment. (See below for more info.)
 
************************************************************************
 
  &&D. WHAT IS THE  c index?
 
  The c index is a general measure of discrimination. Whereas RMS error 
  summarizes how close each pattern's predicted value is to its true 
  outcome, the c index measures the relative goodness of the set of 
  predictions as a whole, by comparing predicted probability of each 
  pattern with that of all other patterns. For this reason, during
  training, a declining c index on an independent test set may be a
  better warning of over-fit (impaired generalization) than RMS error. 
  NevProp offers the choice of c index (BestByCindex set to 1) or RMS 
  error (BestByCindex 0) to determine training set weights that test
  "best" generalize to the testing set.
 
  The computational approach to the c index (ROC curve area) was first 
  reported by Bamber (J Math Psych 1975; 12:387-415), and later
  in the medical literature by Hanley and McNeil (Radiology 1982;
  143:29-36) as the trapezoidal approximation to the area under the ROC
  curve (binary outcomes) and by Harrell for continuous outcomes such
  as duration of survival (JAMA 1982; 247:2543-2546).
 
  The c index is independent of the prior probability of each outcome. 
  Unlike the threshold cutoff for determination of accuracy (the 
  ScoreThreshold parameter in Quickprop), the c index computation
  requires  only that the algorithm produce an ordinally-scaled 
  relative predictive score, not a true probability. 
 
  This version of NevProp generates a c index for only binary outcomes. 
  In this mode, the c index is an estimator of the area under the
  receiver operating characteristic (ROC) curve, which is widely used
  in signal analysis and the evaluation of diagnostic tests. 
  The binary c index estimates the average probability, over all
  pairs of pattern responses with differing true outcome, that the 
  classification algorithm (here, backprop) will assign a higher score 
  to responses with a true binary outcome value = 1 compared to those 
  with a binary outcome value = 0 (or 0.5 and -0.5, respectively). 
 
  For example, if an output node value = 0 means that a patient lived,
  and an output = 1 means death, then the c index estimates the average 
  probability, over all possible pairs of patients with truly differing 
  outcome, that the algorithm will correctly assign a higher predicted 
  mortality probability to those who died than those who survived.
 
************************************************************************
 
  &&E. COMPILING NevProp / Makefile
 
  This code has been successfully compiled under the following:
    o SunOS Release 4.1 (Sun-4/280)
    o HP-UX (HP 9000/750 and 9000/385)
    o Convex OS (Convex C-220)
    o CRAY OS (CRAY-YMP2)
 
    o Macintosh OS (with or without co-processor)
    o DOS (no attempt made to address extended or expanded memory)
 
  To compile the program WITHOUT the c index measure of accuracy, 
  only the "np.c" code is required, compiled under UNIX as:
    "cc np.c -DCALC_CINDEX=0 -lm -o np"
 
  However, if you want to have the c index measure reported and tested
  (the default compilation mode), we suggest you use the Makefile
  provided. All 10 of the following files must then be available:
 
    np.c,    cind.c,   create.c,  rank.c,      create.h,  
    cind.h,  const.h,  do.h,      protoTyp.h,  rank.h
 
  On UNIX workstations, we suggest you compile using the "Makefile" 
  by typing "make". The default Makefile defines rand48() as available. 
  If you don't have rand48(), but do have random(), you can
  define it in the Makefile. Both are better randomization functions
  than the np.c internal default rand(), which is used if neither
  rand48() or random() are defined at compile time. Macintosh and DOS
  executable versions use rand().
   
  If you have GNU C compiler, gcc, on SunOS, you may achieve faster 
  performance compiling with -O2 and -finline, and -static link.
 
 
  NOTE - COMPILING WITH PRE-ANSI C:
 
  NevProp should compile with either an ANSI or traditional C compiler.
  However, unless you compile it with CALC_CINDEX defined as 0, it
  currently requires the equivalent of the ANSI C library function
  qsort(), and the header file stdlib.h.  If your system doesn't have
  one, try creating stdlib.h in the NevProp directory (#include
  "stdlib.h" is used so it need not be placed in a system directory),
  merely containing the declarations below.  You may also need to cast
  the type of the second argument in calls to qsort().
 
  extern void malloc();
  extern void free();
  extern void qsort();
  extern void exit();
 
************************************************************************
 
  &&F  .net FILE and PARAMETER DESCRIPTIONS
 
  The sample 'iris.net' file shows the formatting of a network, or .net
  file. In 'iris.net', you will see 16 operating parameters, followed by
  3 unit number parameters, an output unit range declaration, and 
  'Connectcalls' followed by that number of lines specifying connections
  between ranges of units. These connections can jump across layers.
  The last 2 parameters indicate how many subsequent training followed by 
  testing cases to read in. Here are the parameters, and abbreviated data:

  #
  IBMorDOS 0          UseQuickProp 1        EpochWiseUpdate 1
  
  BestByCindex 0      MinEpochs 200         BeyondBestEpoch 1.5
 
  WeightRange 1.0     HyperErr 0            SigmoidPrimeOffset 0.1
  Epsilon 0.01        SplitEpsilon 1        Momentum 0.1
  WeightDecay -0.001  ScoreThreshold 0.10
  
  QPMaxFactor 1.75    QPModeSwitchThreshold 0.0
  
  Ninputs 4             Nhidden 3           Noutputs 1
  OutputUnitType 2
  Connectcalls 2
  1  4      5  7
  5  7      8  8

  NTrainingPatterns 50   NTestingPatterns 50
 
  5.000 2.000 3.500 1.000    1
  5.900 3.000 4.200 1.500    1
          .
  (46 TrainingPatterns deleted here)
          .
  6.400 2.700 5.300 1.900    0
  6.000 2.700 5.100 1.600    1
  
  5.600 2.900 3.600 1.300    1
  6.800 3.200 5.900 2.300    0
           .
  (46 TestingPatterns deleted here)
           .
  6.300 2.800 5.100 1.500    0
  6.400 3.200 4.500 1.500    1
 

- - - - - - - Some comments on the .net file construction - - - - - - -

o The active part of a network file begins with a pound sign.
  All text before the pound sign is ignored, so you can insert
  information about the data. When creating your own network files, 
  be sure to include that pound sign before the parameters begin 
  (or you will get an error message).

o NevProp assumes that a ".net" extension is appended to the name of
  the network file. Therefore, at the program's prompt, enter only
  the main part of file's name.

o First, you specify the 16 operating parameters (described below).
  We've tried to cluster them in a logical manner - the first
  row are switches for general operation, the second row are
  settings used only when STOPPED TRAINING is invoked interactively
  during NevProp operation, the next 3 rows are setting used in both
  quickprop and regular backprop modes, and the last row's settings
  are used only for quickprop operation.

o Next, you specify the number of input, hidden (total for all
  hidden layers) and output nodes, respectively.

o Next, you declare whether the range of output is -0.5 to +0.5
  (OutputUnitType 1), or 0 to 1 (OutputUnitType 2).

o Next, you declare the number of connection calls you will make on
  subsequent lines, one call for each connection between RANGES of
  nodes. For example, in the iris.net, below, we specify 2
  ConnectCalls. The first call connects each of units 1 through 4
  (the input nodes, describing four features of the flowers) to
  each of units 5 through 7 (we've arbitrarily chosen 3 hidden
  units). The second call connects each of the hidden units to
  the single output node. Notice that if a single output node is
  used, as in this case, it's unit assignment should be LISTED TWICE.

o Next you specify the number of training and testing vectors.
  Here, "NTrainingPatterns 50    NTestingPatterns 50" indicates that
  50 training exemplars will be followed by 50 testing vectors. You
  may use any number of spaces and tabs between vector elements, and
  any number of blank lines between successive vectors.

  NOTE: If you are not providing a set of vectors for testing,
        just enter "NTestingPatterns 0".


- - - - PARAMETERS: the following 6 were created in NevProp: - - - - -
 
IBMorDOS: 1=>keep names of saved files short, for IBM or DOS systems.
          0=>uses extended prefix with "Ninput-Nhidden-Nouput" added. 

UseQuickProp: 1=>Quickprop is preferred weight updating algorithm.
              0=>Do not use QuickProp at all, only gradient descent.

EpochWiseUpdate: 1=>Update the weights after an epoch of 
                    collecting errors.
                 0=>Update weights after each training pattern is
                    presented (at present, this works only for gradient
                    descent - NevProp will force EpochWiseUpdate to 1
                    whenever you set UseQuickProp to 1 in the .net file)

BestByCindex: 1=>Determine best-generalizing weights by highest
                 test set c index calculated for one of the reporting
                 periods.
              0=>Determine best generalization by highest RMS error.
                 [Applies to both quickprop and gradient descent.]  

MinEpochs: When saving best generalizing weights, do not even
           consider stopping before this number of epochs of
           training have occurred.  This is to protect against
           falling into a local minimum of generalization error
           during the early oscillations of network training.
           [Applies to both quickprop and gradient descent.]  

BeyondBestEpoch: When saving best generalizing weights, exceed the epoch
                 at which best generalization occurred by this factor.
                 (E.g., if best generalization occurred at epoch 200 and 
                 BeyondBestEpoch was set to 1.5, NevProp would not stop
                 training (and testing generalization) until epoch 300.
                 This is to protect against prematurely accepting a local
                 minimum of generalization error.)
                 [Applies to both quickprop and gradient descent.]  


- - - PARAMETERS: the following 10 were incorporated from Quickprop - - -
  
WeightRange: +/- range for randomly initialized weights.
             [Applies to both quickprop and gradient descent.]  

HyperErr: 1=> use atanh error function; 0=> use standard sigmoid.
              [Applies to both quickprop and gradient descent.]  

SigmoidPrimeOffset: Add this to sigmoid-prime to kill flat spots.
                    [Applies to both quickprop and gradient descent.]  

Epsilon: Learning rate
         [Applies to both quickprop and gradient descent.]  

SplitEpsilon: 1=> divide Epsilon by number of input nodes (fan-in)
                  before use (i.e., normalize).
              0=> use full Epsilon value for every weight's update.
                 [Applies to both quickprop and gradient descent.]  

Momentum: Momentum term for gradient descent.

WeightDecay: Rate of weight decay - adds this negative portion of previous
             weight to the next set.
            [Applies to both quickprop and gradient descent.]  

ScoreThreshold: (range 0.0 - 1.0) A prediction this close to the true 
                value is considered an accurate bit. It's used only
                in reporting, not training.
               [Applies to both quickprop and gradient descent.]  

QPMaxFactor: In Quickprop, don't jump more than this times last step.

QPModeSwitchThreshold: When running quickprop, it will switch to regular
                       gradient descent inside +/- this threshold.

************************************************************************

&&G. SAMPLE OUTPUT OF RUN on "iris.net"

Below are 2 runs on iris.net: the first run uses the network (named 
'iris_rms.net') with the parameter 'BestByCindex' set to 0 (implies that
RMS error will be used to determine epoch of best test set generalization;
the second, using the same random seed '4444' uses 'iris_c.net' with 
'BestByCindex' set to 1 (implies the c index will be used). 

NOTE: This was run on an HP9000/700 with rand48() available -- your own
      results may differ.

       || Comments are indented and inserted 
       || to the right of double bars, like this.

Here are the results of the first run...

########## Welcome to NevProp, version 1.13, updated May 23,1993 ##########

    Quickprop program by Scott Fahlman (sef+@cs.cmu.edu), translated
    from Lisp to C by Terry Regier.

    Generalization, stopped training, c index, and other enhancements
    by Phil Goodman, David Rosen, and Allen Plummer
    University of Nevada Center for Biomedical Modeling Research

    PLEASE REGISTER to receive update notices - just send a message
    describing your research by Internet (goodman@unr.edu), or by surface
    mail to the address provided in the introduction.

       || You are then prompted for the introductory stuff... 

>>Want the introduction [y=>YES, n=>NO]?[n]

       || Just hit a RETURN to accept the default NO response. 

>>Enter name of network [assumes .net]: iris_rms
  ... Parameters & data were read from "iris_rms.net"

       || The .net file must be in the format described above
       || section &&F, as demonstrated in the sample 'iris.net'. 

>>Enter name of file to copy results [assumes .res]:[(none)] iris_rms
  ...Copy of results will be saved to new "iris_rms.res"

       || Supply a file name prefix if you want to save a copy of the
       || results (.res extension will be added). (if the file already
       || exists the results will be appended). Hit a RETURN if
       || you only want to report results to standard output.

>>Retrieve existing weights?[n]

       || If you had answered "y", you would then be prompted for 
       || the name of an existing .wts file. You will be offered
       || the opportunity to create a prediction file and exit, or 
       || to train starting with the uploaded weights.

  >Use clock time as random seed for weights?[y] n
   >Enter seed for random number generator: 4444

       || Here, we did not accept clock time as the seed, instead 
       || forcing '4444' as the seed for random initial weights. 

>>Retain BEST GENERALIZING weights to stop training?[y]

       || Since we used the default 'y', NevProp will keep track
       || of the weights for the run with the best RMS Error
       || or c index (as specified by "BestByCindex" in the .net file).

       || Answer 'n' if you just want to train the network for
       || a specified number of epochs; you would also be asked
       || whether to run test set results periodically. Only the
       || final weights would be available for saving to a .wts file.

      
>>Enter maximum number of epochs to train: 500

       || You can add on more epochs later, so start low if you are
       || running interactively.  If you have redirected standard
       || input from a text file of responses, used a larger number,
       || because NevProp STOPPED TRAINING routine will control the run.

>>Enter number of epochs between reports:[25]

       || NevProp always offers a default of 20 "peeks" evenly spaced
       || over the maximum number of epochs you selected. Here, we 
       || accepted this default to inspect and judge progress every
       || (500/20) = 25 training epochs.


       || Now it regurgitates information from the .net file...

**********************************************************************
Program started on Sun May 23 16:49:59 1993


Parameters read from file "iris_rms.net"...

TRAINING set patterns: 50
TESTING set patterns:  50

InputUnits: 4 -- HiddenUnits: 3 -- OutputUnits: 1

       || "HiddenUnits" is the sum of all hidden units, regardless
       || of the number of hidden layers you created by your 
       || ConnectCalls in the .net file.

SEED for initial random weights=4444; Using lrand48(),srand48().

       || The seed is reported, in case you want to exactly
       || reproduce a run at a later time.

       || It also tells you which randomizing function is defined.
       || Compile with lrand48() if you have it on your system, 
       || since it is better than random() and the default rand().

IBMorDOS 0              UseQuickProp 1          EpochWiseUpdate 1
BestByCindex 0          MinEpochs 200           BeyondBestEpoch 2
WeightRange 0.1         HyperErr 1              SigmoidPrimeOffset 0.1
Epsilon 0.1             SplitEpsilon 1          Momentum 0.1
Decay -0.001            ScoreThreshold 0.1
MaxFactor 1.75          ModeSwitchThreshold 0

       || You can confirm or later refer to these parameters.

       || Now come the reports; "Epoch 0" is before any training.
       || The TRAINING line refers to progress on the training set,
       || and the TESTING line refers to progress on the test set.
       || If you do not provide a test set (i.e., set NTestingPatterns 
       || to 0 in the .net file), the TESTING line will not appear.

       || If you set "UseQuickProp" to 1 in the .net file, you will
       || find that usually backprop "Did QuickProp on 100.00 %"
       || of patterns in that Epoch, but may be a lower % if the
       || algorithm resorted to regular gradient descent inside the
       || "ModeSwitchThreshold" you set in the .net file.
 
       || If you set UseQuickProp to 0, however, you will always see
       || "Did QuickProp on 0.00 %, GradDesc on 100.00 % of weights."

       || The %correct is based on whether a pattern's prediction
       || is within the ScoreThreshold you specified in the .net file.

       || The RMSError is a 0-1 scaled report of the sum of squared 
       || errors used to train the network in every Epoch. The c 
       || indexes, however, are calculated only at the specified
       || reporting frequency; they are not used in the bp algorithm.

*--------------------------------------------------------------------*
Epoch 0:
 TRAINING:   0.00 %correct ;  RMSErr=0.49907 ;  c index=0.87660
  TESTING:   0.00 %correct ;  RMSErr=0.50001 ;  c index=0.94231
*--------------------------------------------------------------------*
Epoch 25: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:   0.00 %correct ;  RMSErr=0.59240 ;  c index=0.99679
  TESTING:   0.00 %correct ;  RMSErr=0.57189 ;  c index=0.99199
*--------------------------------------------------------------------*
Epoch 50: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:   0.00 %correct ;  RMSErr=0.53762 ;  c index=0.95673
  TESTING:   0.00 %correct ;  RMSErr=0.52335 ;  c index=0.97436
*--------------------------------------------------------------------*
Epoch 75: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:   0.00 %correct ;  RMSErr=0.49574 ;  c index=0.99840
  TESTING:   0.00 %correct ;  RMSErr=0.49066 ;  c index=0.99519
*--------------------------------------------------------------------*
Epoch 100: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:   0.00 %correct ;  RMSErr=0.51052 ;  c index=0.99840
  TESTING:   0.00 %correct ;  RMSErr=0.49574 ;  c index=0.99519
*--------------------------------------------------------------------*
Epoch 125: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  48.00 %correct ;  RMSErr=0.65401 ;  c index=0.99519
  TESTING:  52.00 %correct ;  RMSErr=0.62906 ;  c index=0.99840
*--------------------------------------------------------------------*
Epoch 150: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  36.00 %correct ;  RMSErr=0.23428 ;  c index=0.99359
  TESTING:  34.00 %correct ;  RMSErr=0.27057 ;  c index=0.99038
*--------------------------------------------------------------------*
                                 -
                                 -  (middle reports skipped)
                                 -
*--------------------------------------------------------------------*
Epoch 450: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  92.00 %correct ;  RMSErr=0.13177 ;  c index=0.99679
  TESTING:  86.00 %correct ;  RMSErr=0.23209 ;  c index=0.98878
*--------------------------------------------------------------------*
Epoch 475: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  90.00 %correct ;  RMSErr=0.12971 ;  c index=0.99679
  TESTING:  84.00 %correct ;  RMSErr=0.24482 ;  c index=0.98718

OOO-AHHH! STOPPED TRAINING: Epoch 475 is at least 2 X BestEpoch 225

*--------------------------------------------------------------------*

(Up to) 10 BEST runs on the TEST set, according to RMS error:
  (using training weight SEED=4444)
                                                    RELATIVE
                                                   (1-RMSErr)
                                                   LOW-----HI

 Epoch:   175   RMSErr=0.212893   c index=0.990385 *********
 Epoch:   200   RMSErr=0.218283   c index=0.990385 *******
 Epoch:   225   RMSErr=0.210098   c index=0.990385 **********
 Epoch:   250   RMSErr=0.214798   c index=0.990385 ********
@.................................................................
 Epoch:   300   RMSErr=0.233063   c index=0.987180 ***
 Epoch:   325   RMSErr=0.239659   c index=0.987180 **
 Epoch:   350   RMSErr=0.244285   c index=0.987180 *
 Epoch:   375   RMSErr=0.241557   c index=0.988782 *
 Epoch:   400   RMSErr=0.243707   c index=0.988782 *
@.................................................................
 Epoch:   450   RMSErr=0.232090   c index=0.988782 ****

@..............................=> Skips at least one reporting period.
Training halted after Epoch 475.
*--------------------------------------------------------------------*

       || Notice how, using RMS error minimization as the criterion
       || for best test set generalization, there is an optimum around
       || the 225th epoch, and that the surrounding RMSerrors did
       || fall and rise pretty uniformly (because they are not among the
       || 'top 10', all epochs before 175 had to have larger RMS errors 
       || than epoch 350, the worst displayed in the histogram).

>>Done with training?[y] n

       || We said 'n' in order to force additional training to
       || see if we might have fallen into a local error minimum...

>>Enter number of epochs to FORCE additional training:[500]
>>Enter number of epochs between reports:[25]

       || We accepted the defaults, which will be the same as the
       || most recent prior selections.

***********************************************************************

       || skipping initial re-display of the parameters here; the 
       || reports begin again...

*--------------------------------------------------------------------*
Epoch 500: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  88.00 %correct ;  RMSErr=0.13325 ;  c index=0.99679
  TESTING:  86.00 %correct ;  RMSErr=0.25576 ;  c index=0.98718
*--------------------------------------------------------------------*
Epoch 525: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  90.00 %correct ;  RMSErr=0.12957 ;  c index=0.99679
  TESTING:  84.00 %correct ;  RMSErr=0.24848 ;  c index=0.98718
*--------------------------------------------------------------------*
                                 -
                                 -  (middle reports skipped)
                                 -
*--------------------------------------------------------------------*
Epoch 950: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  90.00 %correct ;  RMSErr=0.12353 ;  c index=0.99840
  TESTING:  86.00 %correct ;  RMSErr=0.25400 ;  c index=0.98718
*--------------------------------------------------------------------*
Epoch 975: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  86.00 %correct ;  RMSErr=0.12999 ;  c index=0.99840
  TESTING:  86.00 %correct ;  RMSErr=0.26993 ;  c index=0.98718

*--------------------------------------------------------------------*

(Up to) 10 BEST runs on the TEST set, according to RMS error:
  (using training weight SEED=4444)
                                                    RELATIVE
                                                   (1-RMSErr)
                                                   LOW-----HI

 Epoch:   175   RMSErr=0.212893   c index=0.990385 *********
 Epoch:   200   RMSErr=0.218283   c index=0.990385 *******
 Epoch:   225   RMSErr=0.210098   c index=0.990385 **********
 Epoch:   250   RMSErr=0.214798   c index=0.990385 ********
@.................................................................
 Epoch:   300   RMSErr=0.233063   c index=0.987180 ***
 Epoch:   325   RMSErr=0.239659   c index=0.987180 *
@.................................................................
 Epoch:   450   RMSErr=0.232090   c index=0.988782 ***
@.................................................................
 Epoch:   625   RMSErr=0.232215   c index=0.987180 ***
@.................................................................
 Epoch:   775   RMSErr=0.233959   c index=0.987180 **
@.................................................................
 Epoch:   925   RMSErr=0.234266   c index=0.987180 **

@..............................=> Skips at least one reporting period.
Training halted after Epoch 975.
*--------------------------------------------------------------------*

       || Notice that we forced 500 more epochs of training, over- 
       || riding the .net parameter that stopped the training initially.

       || The little histogram shows that a few of the later epochs 
       || replace those at the end of the original run, but there was
       || no substantive change around the best generalization peak.
       
>>Done with training?[y]
>>Save BEST test PREDICTIONS to a file?[n] y
   ... Predictions saved to "iris_rms-4-3-1.prd"

>>Save BEST WEIGHTS to a file?[n] y
   ... Weights saved to "iris_rms-4-3-1.wts"

       || If "IBMorDOS" was set to 1 in the .net file, 'iris.prd' and
       || 'iris.wts' would have been created instead, in order to limit
       || the main part of the filename to 8 characters.

       || NevProp will prompt you before over-writing an existing .prd
       || or .wts file.

***********************************************************************
Leaving NevProp...


       || Now, using same random seed, but with the c index as the 
       || best measure (BestByCindex set to 1 in the .net file)...

########## Welcome to NevProp, version 1.13, updated May 23,1993 ##########

       || skipping initial prompts and responses here...

       || Only change is "BestByCindex" now set to 1 in 'iris_c.net'.

**********************************************************************
Program started on Sun May 23 16:48:35 1993


Parameters read from file "iris_c.net"...

TRAINING set patterns: 50
TESTING set patterns:  50

InputUnits: 4 -- HiddenUnits: 3 -- OutputUnits: 1

SEED for initial random weights=4444; Using lrand48(),srand48().

IBMorDOS 0              UseQuickProp 1          EpochWiseUpdate 1
BestByCindex 1          MinEpochs 200           BeyondBestEpoch 2
WeightRange 0.1         HyperErr 1              SigmoidPrimeOffset 0.1
Epsilon 0.1             SplitEpsilon 1          Momentum 0.1
Decay -0.001            ScoreThreshold 0.1
MaxFactor 1.75          ModeSwitchThreshold 0
*--------------------------------------------------------------------*
Epoch 0:
 TRAINING:   0.00 %correct ;  RMSErr=0.49907 ;  c index=0.87660
  TESTING:   0.00 %correct ;  RMSErr=0.50001 ;  c index=0.94231
*--------------------------------------------------------------------*
Epoch 25: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:   0.00 %correct ;  RMSErr=0.59240 ;  c index=0.99679
  TESTING:   0.00 %correct ;  RMSErr=0.57189 ;  c index=0.99199
*--------------------------------------------------------------------*
Epoch 50: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:   0.00 %correct ;  RMSErr=0.53762 ;  c index=0.95673
  TESTING:   0.00 %correct ;  RMSErr=0.52335 ;  c index=0.97436
*--------------------------------------------------------------------*
                                 -
                                 -  (middle reports skipped)
                                 -
*--------------------------------------------------------------------*
Epoch 250: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  84.00 %correct ;  RMSErr=0.13622 ;  c index=0.99359
  TESTING:  82.00 %correct ;  RMSErr=0.21480 ;  c index=0.99038
*--------------------------------------------------------------------*
Epoch 275: Did QuickProp on 100.00 %, GradDesc on   0.00 % of weights.
 TRAINING:  86.00 %correct ;  RMSErr=0.14012 ;  c index=0.99359
  TESTING:  78.00 %correct ;  RMSErr=0.24772 ;  c index=0.98878

OOO-AHHH!! STOPPED TRAINING: Epoch 275 is at least 2 X BestEpoch 125

*--------------------------------------------------------------------*

(Up to) 10 BEST runs on the TEST set, according to c index:
  (using training weight SEED=4444)
                                                    RELATIVE
                                                    C INDEX
                                                   LOW-----HI

 Epoch:    25   RMSErr=0.571894   c index=0.991987 ****
@.................................................................
 Epoch:    75   RMSErr=0.490661   c index=0.995192 *******
 Epoch:   100   RMSErr=0.495743   c index=0.995192 *******
 Epoch:   125   RMSErr=0.629065   c index=0.998397 **********
 Epoch:   150   RMSErr=0.270570   c index=0.990385 **
 Epoch:   175   RMSErr=0.212893   c index=0.990385 **
 Epoch:   200   RMSErr=0.218283   c index=0.990385 **
 Epoch:   225   RMSErr=0.210098   c index=0.990385 **
 Epoch:   250   RMSErr=0.214798   c index=0.990385 **
 Epoch:   275   RMSErr=0.247724   c index=0.988782 *

@..............................=> Skips at least one reporting period.
Training halted after Epoch 275.
*--------------------------------------------------------------------*

>>Done with training?[y]

       || Notice that the c index peaked at Epoch 125, earlier 
       || than RMS Error. This is because the c index is 
       || sensitive to the relative ordering among the predictions,
       || rather than to how calibrated each prediction is with 
       || respect to its individual absolute target.

       || Depending on the structure of your data space, and your
       || classification needs, the c index may give an earlier
       || indication of optimal generalization to a test set.



       || Below, after forcing another 500 training epochs and 
       || skipping the running reports, notice that there is no 
       || change at all among the 'top ten' contenders (note below 
       || that 'Training halted after Epoch 775').

*--------------------------------------------------------------------*

(Up to) 10 BEST runs on the TEST set, according to c index:
  (using training weight SEED=4444)
                                                    RELATIVE
                                                    C INDEX
                                                   LOW-----HI

 Epoch:    25   RMSErr=0.571894   c index=0.991987 ****
@.................................................................
 Epoch:    75   RMSErr=0.490661   c index=0.995192 *******
 Epoch:   100   RMSErr=0.495743   c index=0.995192 *******
 Epoch:   125   RMSErr=0.629065   c index=0.998397 **********
 Epoch:   150   RMSErr=0.270570   c index=0.990385 **
 Epoch:   175   RMSErr=0.212893   c index=0.990385 **
 Epoch:   200   RMSErr=0.218283   c index=0.990385 **
 Epoch:   225   RMSErr=0.210098   c index=0.990385 **
 Epoch:   250   RMSErr=0.214798   c index=0.990385 **
 Epoch:   275   RMSErr=0.247724   c index=0.988782 *

@..............................=> Skips at least one reporting period.
Training halted after Epoch 775.
*--------------------------------------------------------------------*
                                 -
                                 -  (closing prompts skipped)
                                 -
***********************************************************************
Leaving NevProp...


***********************************************************************
&&H. SAMPLE CONTENTS of .wts file

"iris-4-3-1.wts" contains the weights from the best (by RMS error) epoch,
number 225. These weights can be read in at the beginning of another
NevProp session if you (a) wish to begin where you left off, or (b) just 
to make a single predictive pass on a new set of data, without doing any
further training. Here's what the file looks like using a text editor: 

EpochSaved 225
5 0 1.3393
5 1 1.97059
5 2 1.18945
5 3 -2.25401
5 4 -4.48288
6 0 -0.044938
6 1 -0.599122
6 2 -0.243752
6 3 -0.434346
6 4 -0.215725
7 0 -0.567667
7 1 4.7675
7 2 1.92803
7 3 12.1388
7 4 13.5704
8 0 0.727791
8 5 16.8892
8 6 0.00218251
8 7 -3.69998

       || The first line indicates the epoch number at which the
       || subsequent weights were frozen.
 
       || Each line then lists a pair of connected units (hidden 
       || unit index, then the index of each unit connected to it), 
       || followed by its associated weight. Index 0 is for the
       || bias unit. 


***********************************************************************

&&K. SAMPLE CONTENTS of .prd file 
Here's the contents of "iris-4-3-1.prd", 50 predicted test outputs 
based on the saved weights: 

SEQU PRED1 TRUE1
SEQU PRED1 TRUE1
1 0.999992 1.000000
2 0.053856 0.000000
3 0.077359 0.000000
4 0.290028 0.000000
5 0.999994 1.000000
     -
     - (middle lines skipped)
     -
45 0.050382 0.000000
46 0.096132 0.000000
47 0.051701 0.000000
48 0.239320 1.000000
49 0.756762 0.000000
50 0.999877 1.000000

The program first wrote a header row labeling the column variables.
"SEQU" is the sequence of the pattern in the original .net file, 
which can serve as an identifier. "PRED1" indicates that the next column
contains the response of the first output unit. "TRUE1" indicates that 
the last column contains the target outputs of that first unit, as 
provided in the .net file. These space-separated columns should be
readable in free format by most statistical packages.

Notice that Quickprop had problems with patterns 48 and 49; 
these happen to be the same 2 that discriminant analysis misclassifies 
if presented with the same training and testing splits of this data.

Happy net-ing.

* * * * * * * * * * * * END OF "NP_READ.ME" * * * * * * * * * * * * * */
@EOF

chmod 640 NevProp/NP_READ.ME

echo x - NevProp/Makefile
cat >NevProp/Makefile <<'@EOF'
# Makefile for NevProp version 1.03, last modified 13May93
# University of Nevada Center for Biomedical Modeling Research
# goodman@unr.edu

# warning: dependencies on #included .h files are not represented here when
#          they are utilities assumed to be stable.

# Examples:
#   make    (makes the program)
#   make clean   (removes program and object code)
#   make 'PARAMETERS=-DHAVE_RAND48=1' ( machines w/ rand48(); best if avail. )
#   make 'PARAMETERS=-DHAVE_RANDOM=1' ( machines w/ random() but not rand48() )
#         If neither of the above are defined, ran() will be used.
#   make 'PARAMETERS=-DCALC_CINDEX=0' (if c index is NOT desired)
#   make 'OPTIMIZE=-O3'   (for machines that recognize this optimization)
#   make 'OPTIMIZE=-g -DDEBUG'   (debug)

#defaults:
OPTIMIZE = -O
PARAMETERS = -DHAVE_RAND48=1

CFLAGS = $(OPTIMIZE) $(PARAMETERS)
LDLIBS = -lm

TARGET = np

OBJS =	$(TARGET).o \
		cind.o \
		rank.o \
		create.o

$(TARGET): $(OBJS)
		$(CC) $(CFLAGS) -o $(TARGET) $(OBJS) $(LDLIBS)

clean:
	/bin/rm -f $(OBJS) $(TARGET)

@EOF

chmod 640 NevProp/Makefile

echo x - NevProp/np.c
sed 's/^@//' >NevProp/np.c <<'@EOF'
/* * * * * * * * * * * NevProp - - UNIV OF NV CBMR CHANGE LOG * * * * * *
 *
 * 'date'   'vers.' 'who'                'description'
 *
 * 5/26/93   1.15    phg     o Housekeeping; will make first release. 
 *
 * 5/24/93   1.14    phg/dbr o When using c index (defined and Noutputs=1)
 *                             C_INDEX() does the %correct and RMSErr stats
 *                             instead of TRAIN_SUMM() and TEST_SUMMARY().
 *                           o Changed MAKE_PREDICTION() to directly access
 *                             TRAIN() for just the first reporting. Also
 *                             can access it my entering 0 training epochs. 
 *
 * 5/23/93   1.13    phg/dbr o Converted checkEOF from macro to function. 
 *                           o Fixed problem with header row printing
 *                             to .prd file when >1 output unit.
 *                           o Fixed bug in reporting DidGradDesc %.
 *
 * 5/22/93   1.12    phg     o Prints out histogram immediately after every
 *                             training sequence, before prompting for more.
 *                           o Moved "Retain BEST GEN. WTS.." prompt ahead
 *                             of (now) subprompt to run TEST stats w/o
 *                             stopped training (latter unlikely to be used).
 *                           o Supplies default suggestions for ReportFreq
 *                             (initally, 1/20 the epox; on continued 
 *                             training uses previous ReportFreq if legal). 
 *                           o Supplies default suggestion for epox after
 *                             first training sequence (same as previous).
 *                           o Created GET_INTEGER() for 2 defaults above. 
 *                           o Trapped possibility that reportfreq > epox 
 *                             by mistaken entry (allows user to try again).
 *                           o After stopped training, now gives option to
 *                             interactively force more training epochs.
 *                           o Separated out function STOPPED_TRAINING().
 *                           o Separated out function SHUFFLE_10BEST().
 *                           o Separated out function TRAIN_SUMMARY().
 *                           o Gives second try to provide .net & .wts files. 
 *                           o Made fname a global to facilitate second
 *                             tries in filename specification.
 *                             (.net now promtped from GET_NETWORK_CONFIG())
 *                           o Corrected formula for %correct, from Npatterns
 *                             to Npatterns*Noutputs in denominator. 
 *                           o Moved result file prompt to follow .net prompt.
 *
 * 5/21/93   1.11    phg/ajp o Eliminated junk[] using %*s formatting
 *                             in GET_NETWORK_CONFIG & GET_WEIGHTS.
 *                           o Fixed bug in PRINT_10BEST() that 
 *                             permitted epoch 0 error to be considered
 *                             in relative score computation for graphic.
 *
 * 5/20/93   1.10    phg     o Automatically stopped training implemented.
 *                           o New globals added to .net file: MinEpochs
 *                             and BeyondBestEpoch for stop training.
 *                           o Checks for existing x.res, x.wts, and x.prd
 *                             filenames before writing, and offers option
 *                             to assign new names. 
 *                           o Globals rearranged in more logical
 *                             clusters in the .net file.
 *                           o Feeds back names of files read or saved.
 *                           o PRINT_GLOBALS formatting improved.
 *                           o epox=epox+1 (at startup only) to allow for
 *                             Epoch 0 report before Epoch 1 training.
 *                           o Saves and recalls to .wts file the Epoch
 *                             at which it was frozen.
 *                           o Saves completion time to end of res file.
 *
 * 5/18/93   1.04    phg     o Inserted fflush(stdout) after printing
 *                             reports to allow real-time reporting
 *                             to redirected output when not interactively
 *                             connected to terminal. 
 *
 * 5/15/93   1.04    phg     o Re-ordered functions to most common sequence
 *                             called, starting with main(). Corrected
 *                             ordering in Function Declaration list.
 *                           o Declared all functions but main() as static.
 *                           o Changed names of ReadDefault, CIndex, 
 *                             and predict1 to upper case equivalents.
 *                           o Made more uniform the comments format atop
 *                             each function's code.
 *                           o Added int ReportRandSeed to turn off report
 *                             of Random Seed if .wts file was loaded.
 *
 * 5/12/93   1.03    dbr/phg o Macro dbassert installed so that assert() 
 *                             in ACTIVATION() &  ACTIVATION_PRIME() loops
 *                             are called only if DBUG is defined on compile.
 *
 * 5/12/93   1.02    phg     o Incorporated these suggestions by Mike Arras
 *                             (arras@forwiss.uni-erlangen.de) to increase
 *                             speed (he also recommends using gcc with
 *                             -O2 and -finline functions, and -static link):
 *                           o Replaced Errors[] w/ single float 'errors'.
 *                           o ErrorsSums[] initialization removed.
 *                           o In BACKWARD_PASS(), ErrorSums initialization
 *                             starts with FirstHidden.
 *                           o In UPDATE_WEIGHTS(), SplitEpsilon asked only
 *                             once per unit; epsilon set before loop.
 *                           o In FORWARD_PASS() and BACKWARD_PASS(), 
 *                             stop repeated locating of Nconnections[j].
 *                           o Single index referencing in tight loops
 *                             of FORWARD_PASS() and BACKWARD_PASS().
 *
 * 5/11/93   1.01    phg/dbr o In 2 memory allocation lines, changed
 *                             NTrainingPatterns to NTestingPatterns.
 *
 * 5/04/93   1.00    phg     o Local first release. 
 *
 * 5/04/93   0.31    phg     o ACTIVATION() asymptotic thresholds
 *                             removed for more accurate c indexing.
 *
 * 5/04/93   0.31    phg     o #include "ctype.h"
 *
 * 5/03/93   0.30    phg     o Current clock time printed to stdout and
 *                             results file headers.
 *                           o EpochWiseUpdate forced to 1 (w/ message)
 *                             if .net sets it to 0 with UseQuickProp=1.
 *
 * 5/03/93   0.29    phg/dbr o Option added to save stdout results to a
 *                             file, with associated changes in printf's.
 *
 * 5/03/93   0.28    phg     o More housecleaning; inaccessible breaks
 *                             removed from ACTIVATION().
 *
 * 5/02/93   0.27    phg     o Modifications of intro remarks;
 *                             general housekeeping.
 *
 * 5/01/93   0.27    phg     o Moved sRANDOM(seed) and INIT_WEIGHTS() up
 *                             to be called only if random no. requested,
 *                             to allow uploaded old wts to be used.
 *
 * 5/01/93   0.26    ajp     o Updated GET_WEIGHTS for new index format.
 *
 * 4/29/93   0.25    ajp     o Cast math functions (log, exp, fabs and
 *                             sqrt) input as double and outputs as float.
 *                           o Removed limitations in ERRFUN().
 *                           o Finished mods for dynamic memory for PCs.
 *                           o Cast ERRFUN() in TEST_SUMMARY() to void.
 *                           o Made Introduction into a fn.
 *                           o Made a fn for printing globals.
 *                           o Made a fn for printing resulting output.
 *                           o Updated GET_WEIGHTS for new index format.
 *
 * 3/28/93   0.21    phg     o Began memory mods for PCs
 *
 * 3/23/93   0.20    phg     o Added default Unit_type==ASYMSIGMOID
 *                             returns to ACTIVATION_PRIME(value) and
 *                             ACTIVATION(sum) functions.
 *                           o Incremented no. epochs to run by 1
 *                             in order to get final report.
 *                           o exit(1) if .wts file doesn't exist
 *                           o Added option to just save a set of TEST
 *                             predictions to a .prd file; offered if
 *                             ">>Retrieve old weights?" answered yes.
 *                           o Added report of stats w/ above option.
 *
 * 3/22/93   0.19    phg     o Changed plain "cc" default to NO c index,
 *                             but WITH c index if default Makefile is used.
 *                           o Changed getchar() to getchar() so compile is
 *                             OK on Convex machine.
 *
 * 3/20/93   0.18    phg     o Limited no. of "BEST TEST" runs to less than
 *                             10 if fewer reports were actually made.
 *
 * 3/17/93   0.17    dbr     o Use rand48() if HAVE_RAND48 is defined and
 *                             true, else use random() if HAVE_RANDOM is
 *                             defined and true, else use rand() by default.
 *                             Print name of function being used.
 *                           o exit(1) if .net file doesn't exist
 *
 * 3/17/93   0.16    phg     o Abort w/ message if .net file not found.
 *                           o Add more info about c index to intro.
 *                           o Eliminate calls to c index if not defined.
 *
 * 3/16/93   0.15    dbr     o fflush(stdout) before stderr prompts.
 *                           o print error msg again from checkEOF().
 *                           o correct type of argument to time().
 *
 * 3/14/93   0.14    phg     o Added graphic display of ten "best" runs.
 *
 * 3/08/93   0.13    phg     o Added option to save weights (and test set
 *                             predictions) for the "best" generalization as
 *                             judged by either the c index or RMS error.
 *                           o After basic reporting is completed, recalls
 *                             reports of ten "best" generalization runs.
 *
 * 3/08/93           dbr     o Fixed checkEOF() for assert() that is a statement
 *                             rather than expression.  return(0) in tempmain(),
 *                             removed unused variables (suggested by SunOS
 *                             lint *.c).  Declare RANDOM_WEIGHTS() properly at
 *                             top before use (courtesy of HP-UX lint *.c -lm).
 *
 * 3/08/93   0.12    phg     o Added EpochWiseUpdate global switch between
 *                             epoch-wise and pattern-wise weight updating.
 *                             (also added to .net file header listing)
 *                           o Report actual .net filename with regurgitated
 *                             set of global parameters read in.
 *
 * 3/06/93   0.11    phg     o Read user-specified global parameters from the
 *                             .net file header section (must be in fixed
 *                             sequence).
 *                           o Regurgitate parameters read in.
 *                           o Enhanced first line of reporting to specify %
 *                             of patterns using QuickProp vs Grad Descent.
 *                           o Added welcome and quitting remarks.
 *                           o Added UseQuickProp switch to globals read in
 *                             (if==0, forces gradient descent always).
 *                           o Added default (carriage return) response at
 *                             several interactive prompts.
 *
 * 3/05/93   0.10    dbr     o Calculate C Index when there is a single output
 *                             unit. Requires "cind" and "create" utilities
 *                             (and anything they require) unless compiled
 *                             with CALC_CINDEX defined as 0.
 *                           o Use lrand48()
 *                             or random() according to HAVE_RAND48.
 *                           o Check all [f]scanf() calls using checkEOF().
 *                           o Allow user to pre-define MAX_UNITS, MAX_PATTERNS,
 *                             CALC_CINDEX, and HAVE_RAND48
 *                             at compile time without editing this code.
 *                           o Check to make sure we don't exceed MAX_* limits!
 *                             Default MAX_* in this code are chosen to limit
 *                             memory requirements to approx. 1MB.  Unix
 *                             Makefile may substitute greater defaults.
 *                           o Remove unnecessary arrays toutputs, tinputs.
 *                           o Uses assert() from #include "assert.h".
 *
 * 3/04/93   0.08    dbr     o Write prompts to stderr so stdout can be
 *                             redirected separately without prompts.
 *                           o Correct the RMS formula.
 *                           o Print same result format for train, test,
 *                             by using PRINT_RESULTS().
 *                           o Use lrand48 instead of notoriously non-random
 *                             rand()  (change to {,s}random() if {l,s}rand48
 *                             missing).  #include "time.h".
 *                           o Use time() if seed==0 .
 *
 * 3/04/93   0.07    phg/dbr o Report RMS error, etc. on Train set, Test set.
 *
 * 2/21/93   qp5.c   phg     o Inserted interactive option to save predictions.
 *                           o Added a header row to the prediction file, with
 *                             automatic naming of ID, PRED#, and TRUE# vars.
 *                           o Allowed switching to a long, more descriptive
 *                             prediction file name: "filename-I-H-O.prd".
 *                           o Allowed switching to a long, more descriptive
 *                             weight file name: "filename-I-H-O.wts".
 *                           o Added a switch "int IBM" to remove -I-H-O in
 *                             filenames to limit names length for IBM systems.
 *                           o Placed carriage return after each weight saved.
 *
 * 2/21/93   qp4.c   phg     o Changed reporting to specify %accuracy, rather
 *                             than bits wrong.
 *
 * 2/20/93   qp3.c   phg     o Changed SSE reporting to specify RMS error.
 *                           o Inserted global "ints ReportFreq, TestReporting"
 *                             and interactive queries for report frequency
 *                             and whether to TEST with each report.
 *
 * 2/20/93   qp2.c   phg     o Increased #define N 30 to 300
 *                           o Increased float TrainingInputs[200][N] to
 *                                       float TrainingInputs[6000][N].
 *                           o Increased float TrainingOutputs[200][N] to
 *                                       float TrainingOutputs[6000][N].
 *                           o Increased float TestInputs[200][N] to
 *                                       float TestInputs[6000][N].
 *                           o Increased float TestOutputs[200][N] to
 *                                       float TestOutputs[6000][N].
 *
 *   --  Phil Goodman, David Rosen, Allen Plummer
 * * * * * * * * * * * end, UNIV OF NV CBMR CHANGE LOG * * * * * * * * * * */

/* * * * * * * * Quickprop v1 - - Regier & Fahlman CHANGE LOG  * * * * * * *
 * This is Scott Fahlman's quickprop program translated from Common Lisp
 * into C by Terry Regier at the University of California, Berkeley.
 * Netmail address: regier@cogsci.berkeley.edu
 * This version is Quickprop 1 from September, 1988.
 *
 * An example of network setup data is included at the end of this file.
 *
 * The algorithm and some test results are described in Fahlman's paper
 * "Faster-Learning Variations on Back-Propagation: An Empirical Study"
 * in Proceedings of 1988 Connectionist Models Summer School, published
 * by Morgan Kaufmann.
 *
 * Note: the parameter called "mu" in the paper is called "max-factor" here.
 *
 * Changes made to quickprop.c version 1 by N Karunanithi netmail:
 * <karunani@handel.cs.colostate.edu>.
 *
 * Connections can now be specified for multiple ranges of units.
 * For example if you had 3 layers of hidden units and wanted the
 * third layer to have connections to inputs and the second layer,
 * but not the first hidden layer.
 *
 * Bug fix in CONNECT_LAYERS by Richard Dale Romero <rr2p+@andrew.cmu.edu>
 * inserted into the code on September 18, 1991
 *
 * You may specify hidden and output units as sigmoids with ranges
 * of -0.5 to 0.5 (SIGMOIDAL) or from 0.0 to 1.0 (ASYMSIGMOIDAL) in
 * the input file.
 * * * * * * * * * * * END, Regier & Fahlman CHANGE LOG  * * * * * * * */


/* * * * * * * * * * * * NevProp C code begins below * * * * * * * * * */

#include <stdio.h>
#include <math.h>
#include "time.h"
#include "assert.h"
#include "ctype.h" /* declare toupper() */

/* You may pre-define these parameters. Otherwise these defaults apply: -dbr */
#ifndef CALC_CINDEX
#define CALC_CINDEX 1       /* default is to calculate c index */
#endif

#if CALC_CINDEX
#include "create.h"
#include "cind.h"
#endif                    

#ifndef HAVE_RAND48
#define HAVE_RAND48 0       /* Is {l,s}rand48() available? */
#endif

#ifndef HAVE_RANDOM
#define HAVE_RANDOM 0       /* Is {,s}random() available? */
#endif

#if HAVE_RAND48
/* prefer rand48() to random() */
#define lRANDOM lrand48
#define sRANDOM srand48
static char randName[]="lrand48(),srand48()";
#else
#if HAVE_RANDOM
#define lRANDOM random
#define sRANDOM srandom
static char randName[]="random(),srandom()";
#else /* !HAVE_RANDOM */ /* use rand() by default */
void srand();
long lRANDOM()
{ return rand();
}
void sRANDOM(seed)
 unsigned long seed;
{ srand((unsigned int)seed);
}
static char randName[]="rand(),srand()";
#endif /* use rand() by default */
#endif /* HAVE_RAND48 */

/*
 * defines a macro dbassert(expression), same as assert(expression) but checking
 * is *disabled* by default, and enabled only if the symbol DEBUG is defined.
 */
#include "assert.h"
#ifdef DEBUG
#define dbassert(expression) assert(expression)
#else
#define dbassert(expression)
#endif

long lRANDOM();
void sRANDOM();
void *malloc();


/***** Function Declarations *****/
static int   READ_DEFAULT();
static int   GET_INTEGER();
static void  INITIALIZE_GLOBALS();
static void  PRINT_INTRO();
static void  GET_NETWORK_CONFIGURATION();
static void  BUILD_DATA_STRUCTURES();
static void  CONNECT_LAYERS();
static void  GET_WEIGHTS();
static void  MAKE_PREDICTIONS();
static void  INIT_WEIGHTS();
static float RANDOM_WEIGHT();
static void  PRINT_GLOBALS();
static void  TRAIN();
static float C_INDEX();
static void  PREDICT_1();
static void  TRAIN_SUMMARY();
static void  TEST_SUMMARY();
static void  RESHUFFLE_10BEST();
static void  TRAIN_ONE_EPOCH();
static void  FORWARD_PASS();
static float ACTIVATION();
static float ERRFUN();
static int   STOPPED_TRAINING();
static void  CLEAR_SLOPES();
static void  BACKWARD_PASS();
static float ACTIVATION_PRIME();
static void  UPDATE_WEIGHTS();
static void  PRINT_RESULTS();
static void  PRINT_10BEST();
static void  PRINT_GRAPHICS();
static void  TEST_SAVE();
static void  DUMP_WEIGHTS();
static void  FREE_MEMORY();

/***** Constants *****/
#define EXIST "* * * The specified .net file does not exist - aborting * * *\n\n"
#define WTSEXIST "* * * The specified .wts file does not exist - aborting * * *\n\n"
#define BADFREQ "* * * Can't use that report frequency - aborting * * *\n\n"
#define NOPATTERNS "* * * NTraining & NTestingPatterns both 0 in the .net file - aborting * * *\n\n"
#define SIGMOID        1  /* Unit_type is SIGMOID with output=+0.5 to -0.5 */
#define ASYMSIGMOID    2  /* ASYMSIGMOID with output=0.0 to 1.0 */

/*****  Global Variables *****/
char  fname[256];          /* global .net prefix */
int   Epoch;               /* Current epoch number */
float WtRange;             /* Random-init weights in range [-WR,+WR] */
float SigmoidPrimeOffset;  /* Add to sigmoid-prime to kill flat spots */
int   HyperErr;            /* 1=> use atanh error function */
int   SplitEpsilon;        /* 1=> divide epsilon by fan-in before use */
float Epsilon;             /* For grad descent if last step was (almost) 0 */
float Momentum;            /* Normal old momentum term */
float ModeSwitchThreshold; /* Inside thresh, do grad descent; outside, jump. */
float MaxFactor;           /* Don't jump more than this times last step */
float Decay;               /* Weight decay */
int   KeepScore;           /* 1=> accumulate error score for each epoch */
float TotalError;          /* Total output error for one epoch */
float ScoreThreshold;      /* This close to desired value=> bit is correct */
int   TotalErrorBits;      /* Total # bits in epoch that were wrong */
int   DidGradient=0;       /* Total # patterns that did gradient descent */
int   Nunits;              /* Total number of units in net */
int   Ninputs;             /* Number of input units */
int   FirstHidden;         /* Index of 1st hidden unit */
int   Nhidden;             /* Number of hidden units */
int   FirstOutput;         /* Index of 1st output unit */
int   Noutputs;            /* Number of output units */
int   Unit_type;           /* Range of output: 1=> SIGMOID, 2=> ASYMSIGMOID */                        
int   NTrainingPatterns;   /* !! Not in Lisp version.  Needed here. */
int   NTestingPatterns;       /* !! Not in Lisp version.  Needed here. */
FILE  *resfile=NULL;       /* Send results to stdout here also, if desired */
unsigned long seed;
int   ReportFreq=1;     /* #epochs between reports - added 2/20/93 phg */
int   TestReporting;  /* turns on TESTING with reports - added 2/20/93 phg */
int   TrainErrorBits;   /* storage for TRAIN reports */
float TrainSSE;         /* storage for TRAIN reports */
int   TestErrorBits;    /* storage for TESTING reports - added 2/21/93 phg */
float TestSSE;          /* storage for TESTING reports - added 2/21/93 phg */
int   IBM=0;            /* 1=> .prd and .wts filenames shorter, w/o I-H-O */
int   UseQuickProp=1;   /* 0=> forces gradient descent always */
int   WeightCount=0;    /* counts connections, reporting wt update by method */
int   ReportRandSeed=1; /* regurgitate only if weights randomized */
int   EpochWiseUpdate;  /* switches betw per epoch and per pattern update */
int   MinEpochs;        /* min epochs before stopped training by test error */
float BeyondBestEpoch;  /* Stop training when Epoch no. exceeds best Epoch
                              by this multiple */
int   EpochSaved=0;     /* Epoch at which a given .wts file was frozen */ 
int   ForceTraining=0;  /* 1 => override stop-training at end of TRAIN() */
int   PredictOnly=0;    /* 1 => do Epoch 0 only to make predictions only */

float cTest=-1., cTrain=-1.;   /* negative=> unknown */
int   BestEpoch[10];
float BestRMSTest[10];
float BestcTest[10];
float TestRMS=0.0;
int   SaveBestWts=0;
int   BestByCindex=1;
int   SortedEpoch[10];
float SortedRMSTest[10];
float SortedcTest[10];
int   ReportInterv[10]; /* Remembers ReportFreq at time its epoch was run */
int   SortedReportInterv[10];

/* the following pertain to 3/28 mods for Mac & PC dynamic mem allocation */
int   **Connections;    /* C[i][j] lists jth unit projecting to unit i */
float **Weights;        /* W[i][j] holds weight of C[i][j] */
float **DeltaWeights;   /* Change between previous weight and current one */
float **Slopes;         /* Accumulated slope value for each position */
float **PrevSlopes;     /* Similarly, for the last position visited */
float **TrainingInputs; /* Double pointer for array accessing */
float **TrainingOutputs;
float **TestInputs;
float **TestOutputs;
float **BestWeights;    /* Holds weights of best reported test C[i][j] */
float *Outputs;         /* Final output value for each unit */
float *ErrorSums;       /* Total error activation for each unit */
int   *Nconnections;    /* # of INCOMING connections per unit */

#define printboth {printf(printargs); if (resfile) fprintf(resfile, printargs);}

static void checkEOF(scanfReturn) 
     int scanfReturn;
{
  if (scanfReturn==EOF || scanfReturn==0)
  { (void)fprintf(stderr,
		  "ERROR: Premature end-of-file encountered on input.\n");
    (void)exit(1);
  }
}


main()
{
 int j;
 int epox, epoxtrue, response, defaultfreq;
 char rname[256], realresfname[260], forcedepochs[10], reportepochs[10];
 char *starttime, *completedtime;
 time_t now;
 struct tm *date;
 
 now = time( NULL );
 date = localtime( &now );
 starttime = asctime( date );

 for (j=0; j<10; j++) forcedepochs[j] = '0';
 for (j=0; j<10; j++) reportepochs[j] = '0';

 epox=34;  /* arbitrary positive */

 INITIALIZE_GLOBALS();

 fflush(stdout);
 fprintf(stderr, "\n\n");
 fprintf(stderr, "\
########## Welcome to NevProp, version 1.15, updated May 26,1993 ##########\n\n");

 fprintf(stderr, "\
    Quickprop program by Scott Fahlman (sef+@cs.cmu.edu), translated\n\
    from Lisp to C by Terry Regier.\n\
\n\
    Generalization, stopped training, c index, and other enhancements\n\
    by Phil Goodman, David Rosen, and Allen Plummer\n\
    University of Nevada Center for Biomedical Modeling Research\n\
\n\
    PLEASE REGISTER to receive update notices - just send a message\n\
    describing your research by Internet (goodman@unr.edu), or by surface\n\
    mail to the address provided in the introduction.\n\n");

 fprintf(stderr, ">>Want the introduction [y=>YES, n=>NO]?[n] ");
 response=READ_DEFAULT(0);
 if(response) PRINT_INTRO();

 GET_NETWORK_CONFIGURATION();

 if (NTrainingPatterns<1 && NTestingPatterns <1)
  {fprintf(stderr, NOPATTERNS);  exit(1);}

 /* Get result file name */
 fflush(stdout);
 fprintf(stderr, "\n>>Enter name of file to copy results [assumes .res]:[(none)] ");
 gets(rname);
 if (rname[0] != '\0') 
 {
  sprintf(realresfname, "%s.res", rname);
  
  resfile=fopen(realresfname, "r");
  
  if (!resfile)
  {
   fflush(stdout);
   fprintf(stderr, "  ...Copy of results will be saved to new \"%s\"\n", realresfname); 
   resfile=fopen(realresfname, "a");
   assert(resfile);
  }  
  else
  {
   fclose(resfile);
   resfile=fopen(realresfname, "a");
 
   fflush(stdout);
   fprintf(stderr, "  ...Copy of NevProp results will be appended to existing \"%s\"\n", realresfname); 
  }
  fprintf(resfile, "\n\n######################### Welcome to NevProp #########################\n");
  fprintf(resfile, "NevProp started on %s", starttime);
 }
 
 fflush(stdout);
 fprintf(stderr,"\n>>Retrieve existing weights?[n] ");
 
 if (READ_DEFAULT(0)) /* yes, retrieve them... */
 {
  ReportRandSeed=0;
  GET_WEIGHTS();

  if (NTestingPatterns > 0)  
  {
   fflush(stdout);
   fprintf(stderr, "\n  >Create a .prd file from the TEST set using these wts?[n] ");

   if (READ_DEFAULT(0))  MAKE_PREDICTIONS(); /* which will exit the program */
  }
 }
 else  /* if not using old weights, start up the random number generator */
 {
  fflush(stdout);
  fprintf(stderr, "  >Use clock time as random seed for weights?[y] ");
  
  if (!READ_DEFAULT(1))
  {
   fprintf(stderr, "   >Enter seed for random number generator: ");
   checkEOF(scanf("%lu", &seed));
   while (getchar() !='\n');
  }
  sRANDOM(seed);
  INIT_WEIGHTS();
 }

 if (NTestingPatterns > 0 && NTrainingPatterns>0)
 {
  fflush(stdout);
  fprintf(stderr, "\n>>Retain BEST GENERALIZING weights to stop training?[y] ");
  if (READ_DEFAULT(1))
  {
   SaveBestWts=1;
   TestReporting=1;

   for (j=0; j<10; j++)   /* Initialize best value array */
   {
    BestEpoch[j]=-1;
    BestRMSTest[j]=-1.;
    BestcTest[j]=-1.;
    ReportInterv[j]=-1;
   }
  }
  else 
  {
   SaveBestWts=0;

   fflush(stdout);
   fprintf(stderr, "  >Run TEST set with each report (but no stopped training)?[y] ");
   if (READ_DEFAULT(1))  TestReporting=1;
   else TestReporting=0;
  }
 }

 /* Train the sucker. */
 while (epox !=0 && NTrainingPatterns>0)
 {
  if (Epoch==0 || Epoch==EpochSaved)
  {
    fflush(stdout);
    fprintf(stderr, "\n>>Enter maximum number of epochs to train: ");
    checkEOF(scanf("%d", &epox));

    if (epox==0 && NTestingPatterns>0)
    {
     while (getchar() !='\n');
     MAKE_PREDICTIONS(); /* which will then exit the program */
    }

    epoxtrue=epox; /* used for error check and default epox if continuing */

    epox=epox + 1;
    while (getchar()!='\n') ;
  }
  else
  {
   for (j=0; j<10; j++)   /* (Re)initialize sorted value array */
   {
    SortedEpoch[j]=-1;
    SortedRMSTest[j]=-1.;
    SortedcTest[j]=-1.;
    SortedReportInterv[j]=-1;
   }

   PRINT_10BEST(); /* display the little histogram at this point */

   fflush(stdout);
   fprintf(stderr,"\n>>Done with training?[y] ");
   if (READ_DEFAULT(1)) epox=0;
   else
   {
    if (SaveBestWts==1)
    {
     ForceTraining=1; /* overrides the stop-training at end of TRAIN() */
     fflush(stdout); 
     fprintf(stderr, ">>Enter number of epochs to FORCE additional training:[%d] ", epoxtrue);
     if (GET_INTEGER(forcedepochs) == 1) epox = atoi(forcedepochs);
     else epox = epoxtrue;
    }
    else  /* just train on thru to completion of the number of epochs desired */
    {
     fflush(stdout); 
     fprintf(stderr, ">>Enter number of epochs to train:[%d] ", epoxtrue);
     if (GET_INTEGER(forcedepochs) == 1) epox = atoi(forcedepochs);
     else epox = epoxtrue;
    }
    epoxtrue=epox; /* used for error check and default epox if continuing */
    epox=epox + 1; /* allows the initial (Epoch 0 or EpochSaved) to report out */
   }
  }

  if (epox !=0)   /* we need the ReportFreq */
  {
  
   /* if the first time thru, or low new epox, suggest 20 even steps */
   if (Epoch==0 || Epoch==EpochSaved || ReportFreq>epoxtrue) 
   {
    defaultfreq=(int)(0.05*epoxtrue);
    if (defaultfreq==0) defaultfreq=1;
   }
   else  defaultfreq=ReportFreq; /* otherwise, suggest previous ReportFreq step */
   
   fflush(stdout);
   fprintf(stderr, ">>Enter number of epochs between reports:[%d] ", defaultfreq);
   if (GET_INTEGER(reportepochs) == 1) ReportFreq = atoi(reportepochs);
   else ReportFreq = defaultfreq;
   
   if ((ReportFreq>epoxtrue) || ReportFreq==0) /* try again if illegal */
   {
    fprintf(stderr, "\n\
!! Epochs BETWEEN reports must be < total epochs and non-zero; try again...");
   fflush(stdout);
   fprintf(stderr, "\n>>Enter number of epochs BETWEEN reports: ");
   checkEOF(scanf("%d", &ReportFreq));
   while (getchar()!='\n') ;

   if ((ReportFreq>epoxtrue) || ReportFreq==0) {fprintf(stderr, BADFREQ);  exit(1);};
  }

   if (Epoch==0)  fprintf(stderr, "\n\
**********************************************************************\n\
Program started on %s\n", starttime);

   PRINT_GLOBALS(stdout, fname);  /* Regurgitate global variables */
   fflush(stdout);
   if (resfile)  PRINT_GLOBALS(resfile, fname);
   TRAIN(epox);
  }
 }

 /* Test the sucker. */
 if (SaveBestWts)
 {
  fflush(stdout);
  fprintf(stderr,">>Save BEST test PREDICTIONS to a file?[n] ");
  response= READ_DEFAULT(0);
 }
 else if (NTestingPatterns>0)
 {
  fflush(stdout);
  fprintf(stderr,">>Save FINAL test PREDICTIONS to a file?[n] ");
  response= READ_DEFAULT(0);
 }
 if (response)  TEST_SAVE();

 if (SaveBestWts)
 {
  fflush(stdout);
  fprintf(stderr,">>Save BEST WEIGHTS to a file?[n] ");
  response=READ_DEFAULT(0);
 }
 else if (NTrainingPatterns>0 || NTestingPatterns>0)
 {
  fflush(stdout);
  fprintf(stderr,">>Save FINAL WEIGHTS to a file?[n] ");
  response= READ_DEFAULT(0);
 }
 if (response) DUMP_WEIGHTS ();

 fflush(stdout);

 if (resfile)
 {
  now = time( NULL );
  date = localtime( &now );
  completedtime = asctime( date );

  fprintf(resfile, "\
**********************************************************************\n\
NevProp completed on %s", completedtime);
 }
 fclose(resfile);
 
 fflush(stdout);
 fprintf(stderr,"\n\
**********************************************************************\n\
Leaving NevProp...\n\n");
 
 FREE_MEMORY();
 return(0);
}
/****** END OF main() ******/


/*
 *  READ_DEFAULT():
 *  Return the appropriate response to a carriage return at a prompt.
 */
static int READ_DEFAULT(Defresponse)
 int Defresponse;
{
 char TempString[256];
 gets(TempString);

 switch(toupper(TempString[0]))
 {
  case 'Y':
  return(1);

  case 'N':
  return(0);

  case '\n':
  default:
  return(Defresponse);
 }
}
/****** END OF READ_DEFAULT() ******/


/*
 *  GET_INTEGER():
 *  Returns an array (up to 10 digits) and provides
 *  a default carriage return switch.
 */
static int GET_INTEGER(line)
 char *line;
{
 char cc;

 if ((cc = getchar( )) !='\n')
 {  
  *line = cc;
  line++;

  while (( cc = getchar( )) !='\n')
  {
   *line = cc;
   line++;
  }
  *line = 0;
  return(1);
 }
 else
 {
  return(2);
 }
}


/*
 *  INITIALIZE_GLOBALS():
 *  Do just that.
 */
static void INITIALIZE_GLOBALS()
{
 Epoch=0;
 KeepScore=0.0;  /* intermittenlty gets reset to 1 for reporting  */
 TotalError=0.0;
 ScoreThreshold=0.35;
 TotalErrorBits=0;
 TestReporting=1;
 seed=time((void *)0); /* set the seed to the current time for default setting */
}
/****** END OF INITIALIZE_GLOBALS() ******/


/*
 *  PRINT_INTRO():
 *  Display introductory remarks if desired.
 */
static void PRINT_INTRO()
{
 fflush(stdout);
 fprintf(stderr,"\
\n\
 NevProp is a general purpose back-propagation program written in C for\n\
 UNIX, Macintosh, and DOS. The original version was Quickprop 1.0 by\n\
 Scott Fahlman, as translated from Common Lisp into C by Terry Regier.\n\
\n\
 The quickprop algorithm itself was not substantively changed, but we\n\
 inserted options to force gradient descent (per-epoch or per-pattern)\n\
 and added generalization monitoring (to stop training), the c index,\n\
 and interface enhancements.\n\n");

 fprintf(stderr, " Press RETURN to continue.\n\n");
 while (getchar()!='\n') ;

 fprintf(stderr, "\
* * * * * * * * * * * * * * CURRENT FEATURES * * * * * * * * * * * * * *\n\n\
o UNLIMITED number of input training and testing PATTERNS;\n\n\
o UNLIMITED number of input, hidden, and output UNITS;\n\n\
o Arbitrary CONNECTIONS among the various layers' units;\n\n\
o Clock-time or user-specified RANDOM SEED for initial random weights;\n\n\
o Choice of regular GRADIENT DESCENT or QUICKPROP;\n\n\
o Choice of LOGISTIC or TANH activation functions;\n\n\
o Choice of PER-EPOCH or PER-PATTERN (stochastic) weight updating;\n\n");

 fprintf(stderr, " Press RETURN to continue.\n\n");
 while (getchar() !='\n');

 fprintf(stderr, "\
o GENERALIZATION to a test dataset;\n\n\
o AUTOMATICALLY STOPPED TRAINING based on generalization;\n\n\
o RETENTION of best-generalizing weights and predictions;\n\n\
o Simple but useful bar GRAPH to show smoothness of generalization;\n\n\
o SAVING of results to a file while working interactively;\n\n\
o SAVING of weights file and reloading for continued training;\n\n\
o PREDICTION-only on datasets by applying an existing weights file;\n\n");

 fprintf(stderr, "\
o In addition to RMS error, the concordance, or c index is displayed.\n\
  The c index shows the correctness of the RELATIVE ordering\n\
  of predictions AMONG the cases; ie, it considers all possible PAIRS\n\
  of vectors. This statistic is identical to the area under the\n\
  receiver operating characteristic (ROC) curve, widely used in\n\
  technology assessment. (See the NP_READ.ME for more info.)\n\n");

 fprintf(stderr, " Press RETURN to continue.\n\n");
 while (getchar()!='\n');

 fprintf(stderr, "\
 NevProp is distributed as freeware with the understanding that\n\
 the code will not be incorporated into commercial applications.\n\
 No warrantee of performance or system compatibility is implied.\n\n");

 fprintf(stderr, "\
 Questions and comments regarding the quickprop algorithm itself should\n\
 be directed to Scott Fahlman at sef+@cs.cmu.edu. Issues pertaining\n\
 to NevProp code or functionality should come our way.\n\n\
 Phil Goodman, David Rosen, and Allen Plummer\n\
 University of Nevada Center for Biomedical Modeling Research\n\
 Internet:  goodman@unr.edu\n\
 Surface :  Washoe Medical Center, Room H-166, 77 Pringle Way,\n\
            Reno, NV  89520   702-328-4867    FAX:702-328-4111\n\n");

 fprintf(stderr, " Press RETURN to continue.\n\n");
 while (getchar()!='\n');

 fprintf(stderr, "\
 The latest version of NevProp is available in C by anonymous ftp from\n\
 \"unssun.scs.unr.edu\" in the directory \"pub/goodman/nevpropdir\".\n\n");
 
 fprintf(stderr, "\
 Macintosh and DOS executable versions may be available in the unr\n\
 anonymous ftp above. If you are unable to transfer the files, just\n\
 send a formatted 3.5 inch floppy disk (specifying Mac or DOS) with\n\
 a self-addressed stamped return envelope to Phil Goodman at the address\n\
 above, including a brief description of yourself and your research.\n\n");

 fprintf(stderr, "\
 Before reporting bugs, please check the unr machine for the \"nevprop.bug\"\n\
 file for bug reports, fixes, and upgrading plans.\n\n");
}
/****** END OF  PRINT_INTRO() ******/


/*
 *  GET_NETWORK_CONFIGURATION():
 *  Read in parameters and allocate memory for data, initialize the network.
 */
static void GET_NETWORK_CONFIGURATION()
{
 FILE  *infile, *fopen();
 char realfname[80];
 char c;
 int  temp[10], i, j, connect;
 float tempfloat, *floatPtr;

 fflush(stdout);
 fprintf(stderr, "\n>>Enter name of network [assumes .net]: ");
 checkEOF(scanf("%s", fname));
 while (getchar() !='\n');

 sprintf(realfname, "%s.net", fname);
 infile=fopen(realfname, "r"); 
 if (!infile)
 {
  fflush(stdout);
  fprintf(stderr,"\n!!That .net file wasn't found; try again...\n");

  fflush(stdout);
  fprintf(stderr, "\n>>Enter name of network [assumes .net]: ");
  checkEOF(scanf("%s", fname));
  while (getchar() !='\n');
 
  sprintf(realfname, "%s.net", fname);
  infile=fopen(realfname, "r");

  if (!infile) {fprintf(stderr, EXIST);  exit(1);}
 }

 c='c';           /* Discard leading comments */
 while (c !='#')  checkEOF(fscanf(infile, "%c", &c));

 /* Get global variables */
 checkEOF(fscanf(infile, "%*s %d %*s %d %*s %d",
  &IBM, &UseQuickProp, &EpochWiseUpdate));

 if (UseQuickProp==1 && EpochWiseUpdate==0) 
 {
  EpochWiseUpdate=1;
  fprintf(stderr,"\n\nNOTE: EpochWiseUpdate was set to 0 in \"%s.net\".\n\
      However, EpochWiseUpdate must be used if UseQuickProp is requested;\n\
      therefore, EpochWiseUpdate is reset to 1 for the run below.\n\n", fname);
 }

 checkEOF(fscanf(infile, "%*s %d %*s %d %*s %f",
  &BestByCindex, &MinEpochs, &BeyondBestEpoch));

 checkEOF(fscanf(infile, "%*s %f %*s %d %*s %f",
  &WtRange, &HyperErr, &SigmoidPrimeOffset));
 checkEOF(fscanf(infile, "%*s %f %*s %d %*s %f",
  &Epsilon, &SplitEpsilon, &Momentum));
 checkEOF(fscanf(infile, "%*s %f %*s %f",
  &Decay, &ScoreThreshold));
 checkEOF(fscanf(infile, "%*s %f %*s %f",
  &MaxFactor, &ModeSwitchThreshold));

 /* Get numbers of inputs, hidden units, and output units */
 checkEOF(fscanf(infile, "%*s %d %*s %d %*s %d",
  &temp[0], &temp[1], &temp[2]));
 BUILD_DATA_STRUCTURES(temp[0], temp[1], temp[2]);

 /* Get the type units used in hidden and outpt layers. */
 checkEOF(fscanf(infile, "%*s %d ", &temp[0]));
 if (temp[0]==1)  Unit_type=SIGMOID;
 else if (temp[0]==2)  Unit_type=ASYMSIGMOID;

 /* Connect layers. */
 checkEOF(fscanf(infile, "%*s %d", &connect));

 for (i=0; i<connect; i++)      /* Reading CONNECT_LAYERS lines */
 {
  checkEOF(fscanf(infile, "%d %d %d %d",
   &temp[0], &temp[1], &temp[2], &temp[3]));
  CONNECT_LAYERS(temp[0], temp[1], temp[2], temp[3]);
 }

 /* Read in number of training patterns */
 checkEOF(fscanf(infile, "%*s %d %*s %d",
  &NTrainingPatterns, &NTestingPatterns)); /* mod 3/28 */

 /* Now get memory for training and testing patterns */
 TrainingInputs=(float **)malloc(((NTrainingPatterns +1) * sizeof(floatPtr)));
 assert(TrainingInputs);
 for(i=0;i<NTrainingPatterns;i++)
 {
  TrainingInputs[i]=(float *)malloc(((Ninputs+1) * sizeof(float)));
  assert(TrainingInputs[i]);
 }

 TrainingOutputs=(float **)malloc(((NTrainingPatterns+1) * sizeof(floatPtr)));
 assert(TrainingOutputs);
 for(i=0;i<NTrainingPatterns;i++)
 {
  TrainingOutputs[i]=(float *)malloc(((Noutputs+1) * sizeof(float)));
  assert(TrainingOutputs[i]);
 }

 TestInputs=(float **)malloc(((NTestingPatterns+1) * sizeof(floatPtr)));
 assert(TestInputs);
 for(i=0;i<NTestingPatterns;i++)
 {
  TestInputs[i]=(float *)malloc(((Ninputs+1) * sizeof(float)));
  assert(TestInputs[i]);
 }

 TestOutputs=(float **)malloc(((NTestingPatterns+1) * sizeof(floatPtr)));
 assert(TestOutputs);
 for(i=0;i<NTestingPatterns;i++)
 {
  TestOutputs[i]=(float *)malloc(((Noutputs+1) * sizeof(float)));
  assert(TestOutputs[i]);
 }

 /* Read patterns into their respective memory blocks */
 for (i=0; i<NTrainingPatterns; i++)
 {
  for (j=0; j<Ninputs; j++)
  {
   checkEOF(fscanf(infile, "%f", &tempfloat));
   TrainingInputs[i][j]=tempfloat;
  }
  for (j=0; j<Noutputs; j++)
  {
   checkEOF(fscanf(infile, "%f", &tempfloat));
   TrainingOutputs[i][j]=tempfloat;
  }
 }

 /* Read in the test patterns */
 for (i=0; i<NTestingPatterns; i++)
 {
  for (j=0; j<Ninputs; j++)
  {
   checkEOF(fscanf(infile, "%f", &tempfloat));
   TestInputs[i][j]=tempfloat;
  }
  for (j=0; j<Noutputs; j++)
  {
   checkEOF(fscanf(infile, "%f", &tempfloat));
   TestOutputs[i][j]=tempfloat;
  }
 }
 fclose(infile);
 
 fflush(stdout);
 fprintf(stderr, "  ... Parameters & data were read from \"%s\"\n", realfname); 
}
/****** END OF GET_NETWORK_CONFIGURATION() ******/


/*
 *  BUILD_DATA_STRUCTURES():
 *  Set up the network connections and allocate necessary memory.
 */
static void BUILD_DATA_STRUCTURES(ninputs, nhidden, noutputs)
 int ninputs, nhidden, noutputs;
{
 int i, *intPtr;
 float *floatPtr;

 Nunits     =1 + ninputs + nhidden + noutputs;
 Ninputs    =ninputs;
 FirstHidden=1 + ninputs;
 Nhidden    =nhidden;
 FirstOutput=1 + ninputs + nhidden;
 Noutputs   =noutputs;

 Connections=(int **)malloc((Nunits * sizeof(intPtr)));
 assert(Connections);
 for(i=0;i<Nunits;i++)
 {
  Connections[i]=(int *)malloc((Nunits * sizeof(int)));
  assert(Connections[i]);
 }

 Weights=(float **)malloc((Nunits * sizeof(floatPtr)));
 assert(Weights);
 for(i=0;i<Nunits;i++)
 {Weights[i]=(float *)malloc((Nunits * sizeof(float)));
 assert(Weights[i]);}

 BestWeights=(float **)malloc((Nunits * sizeof(floatPtr)));
 assert(BestWeights);
 for(i=0;i<Nunits;i++)
 {BestWeights[i]=(float *)malloc((Nunits * sizeof(float)));
 assert(BestWeights[i]);}

 DeltaWeights=(float **)malloc((Nunits * sizeof(floatPtr)));
 assert(DeltaWeights);
 for(i=0;i<Nunits;i++)
 {DeltaWeights[i]=(float *)malloc((Nunits * sizeof(float)));
 assert(DeltaWeights[i]);}

 Slopes=(float **)malloc((Nunits * sizeof(floatPtr)));
 assert(Slopes);
 for(i=0;i<Nunits;i++)
 {Slopes[i]=(float *)malloc((Nunits * sizeof(float)));
 assert(Slopes[i]);}

 PrevSlopes=(float **)malloc((Nunits * sizeof(floatPtr)));
 assert(PrevSlopes);
 for(i=0;i<Nunits;i++)
 {PrevSlopes[i]=(float *)malloc((Nunits * sizeof(float)));
 assert(PrevSlopes[i]);}

 Outputs=(float *)malloc((Nunits * sizeof(float)));
 assert(Outputs);

 ErrorSums=(float *)malloc((Nunits * sizeof(float)));
 assert(ErrorSums);

 Nconnections=(int *)malloc((Nunits * sizeof(int)));
 assert(Nconnections);

 for (i=0; i<Nunits; i++)    Outputs[i]=0.0;
 for (i=0; i<Nunits; i++)    Nconnections[i]=0;

 Outputs[0]=1.0;        /* The bias unit */
}
/****** END OF BUILD_DATA_STRUCTURES() ******/


/*
 *  CONNECT_LAYERS():
 *  Build a connection from every unit in range1 to every unit in range2.
 *  Also add a connection from the bias unit (unit 0) to every unit in
 *  range2.  Set up random weights on links.
 */
static void CONNECT_LAYERS(start1, end1, start2, end2)
 int start1, end1, start2, end2;
{
 int i, j, k;

 for (i=start2; i<=end2; i++)
 {
  if(Nconnections[i]==0)
  {
   Nconnections[i]  +=1;
   Connections[i][0]=0;
   DeltaWeights[i][0]=0.0;
   Slopes[i][0]=0.0;
   PrevSlopes[i][0]=0.0;
   k=1;
  }
  else  k=Nconnections[i];

  for (j=start1; j<=end1; j++)
  {
   Nconnections[i]  +=1;
   Connections[i][k]=j;
   DeltaWeights[i][k]=0.0;
   Slopes[i][k]=0.0;
   PrevSlopes[i][k]=0.0;
   k++;
  }
 }
}
/****** END OF CONNECT_LAYERS() ******/


/*
 *  GET_WEIGHTS():
 *  Retrieve existing weights from a .wts file if requested.
 */
static void GET_WEIGHTS()
{
 FILE  *infile, *fopen();
 int  i, j, k, l;
 float inweight;
 char fwtname[256], tempname[270], realfname[280];

 /* generate full default .wts file name */
 if (IBM==0)
  sprintf(tempname, "%s-%d-%d-%d", fname, Ninputs, Nhidden, Noutputs);
 else
  sprintf(tempname, "%s", fname);

 fprintf(stderr,"  >Enter name of weights file [assumes .wts]:[%s] ", tempname);
 gets(fwtname);
 if (fwtname[0] != 0)  sprintf(realfname, "%s.wts", fwtname);
 else sprintf(realfname, "%s.wts", tempname);

 infile=fopen (realfname, "r");
 if (!infile)
 {
  fflush(stdout);
  fprintf(stderr,"\n!!That .wts file wasn't found; try again...\n");

  fprintf(stderr,"  >Enter name of weights file [assumes .wts]:[%s] ", tempname);
  gets(fwtname);
  if (fwtname[0] != 0)  sprintf(realfname, "%s.wts", fwtname);
  else sprintf(realfname, "%s.wts", tempname);

  infile=fopen (realfname, "r");

 
  if (!infile) {fprintf(stderr, WTSEXIST);  exit(1);}
 }
 /* Read in Epoch at which the wts were frozen */
 checkEOF(fscanf(infile, "%*s %d", &EpochSaved));

 for (i=0; i<Nunits; i++)
  for (j=0; j<Nunits; j++)  Weights[i][j]=0.0;    /* Zeroed weights */

 /* Read in Weights based on .net file construction */
 for (j=FirstHidden; j<Nunits; j++)
  for (i=0; i<Nconnections[j]; i++)
  {
   checkEOF(fscanf(infile, "%d %d %f", &k, &l, &inweight));
   /* Check array elements for proper placement of weight values */
   if ((j==k) && (l==Connections[j][i]))  Weights[j][i]=inweight;
   assert(l==Connections[j][i]);
  }
 
 Epoch=EpochSaved;

 fclose (infile);

#define printargs\
 "... Weights saved at Epoch %d retrieved from \"%s\"\n", EpochSaved, realfname
 printboth;
#undef printargs

}
/****** END OF  GET_WEIGHTS() ******/


/*
 *  MAKE_PREDICTIONS():
 *  This can be called from the prompt offered after uploading existing
 *  .wts file, or by selecting 0 epochs at that prompt (in which case
 *  predictions based on initial random wts can be made).
 */
static void MAKE_PREDICTIONS()
{
 PredictOnly=1; /* for report formatting */
 TRAIN(1);  /* so no training, just reporting happens */
 TEST_SAVE();  /* now save the predictions */

 fflush(stdout);
 fprintf(stderr,"\n\
**********************************************************************\n\
Leaving NevProp...\n\n");

 FREE_MEMORY();

 /* ...and terminate the program */
 exit(0);
}
/****** END OF MAKE_PREDICTIONS() ******/


/*
 *  INIT_WEIGHTS():
 *  For each connection, select a random initial weight between WtRange
 *  and its negative.  Clear delta and previous delta values.
 */
static void INIT_WEIGHTS()
{
 int i, j;

 for (i=0; i<Nunits; i++)
 {
  for (j=0; j<Nconnections[i]; j++)
  {
   Weights[i][j]=RANDOM_WEIGHT(WtRange);
   DeltaWeights[i][j]=0.0;
   Slopes[i][j]=0.0;
   PrevSlopes[i][j]=0.0;
  }
 }
}
/****** END OF INIT_WEIGHTS() ******/


/*
 *  RANDOM_WEIGHT():
 *  Return a float between -range and +range.
 */
static float RANDOM_WEIGHT(range)
 float range;
{
 return ((float)(range*(lRANDOM()%1000/500.0)) - range);
}
/****** END OF RANDOM_WEIGHT() ******/


/*
 *  PRINT_GLOBALS():
 *  Regurgitate parameters set from the .net file.
 */
static void PRINT_GLOBALS(stream, filename)
 char *filename;
 FILE *stream;
{
 fflush(stdout);
 fprintf(stream, "\nParameters read from file \"%s.net\"...\n\n",filename);

 fprintf(stream, "TRAINING set patterns: %d\n", NTrainingPatterns);

 fprintf(stream, "TESTING set patterns:  %d\n\n", NTestingPatterns);

 fprintf(stream, "InputUnits: %d -- HiddenUnits: %d -- OutputUnits: %d\n\n",
          Ninputs, Nhidden, Noutputs);

 if (ReportRandSeed==1)
 fprintf(stream, "SEED for initial random weights=%lu; Using %s.\n\n", seed, randName);

 fprintf(stream, "IBMorDOS %d\t\tUseQuickProp %d\t\tEpochWiseUpdate %d\n",
                  IBM, UseQuickProp, EpochWiseUpdate);

 fprintf(stream, "BestByCindex %d\t\tMinEpochs %d\t\tBeyondBestEpoch %G\n",
                  BestByCindex, MinEpochs, BeyondBestEpoch);

 fprintf(stream, "WtRange %G\t\tHyperErr %d\t\tSigmoidPrimeOffset %G\n",
                  WtRange, HyperErr, SigmoidPrimeOffset);

 fprintf(stream, "Epsilon %G\t\tSplitEpsilon %d\t\tMomentum %G\n",
                  Epsilon, SplitEpsilon, Momentum);

 fprintf(stream, "Decay %G\t\tScoreThreshold %G\n",
                  Decay, ScoreThreshold);


 fprintf(stream, "MaxFactor %G\t\tModeSwitchThreshold %G\n",
                  MaxFactor, ModeSwitchThreshold);
}
/****** END OF PRINT_GLOBALS() ******/


/*
 *  TRAIN():
 *  Train the network for the specified number of epochs, printing out
 *  performance stats every 'ReportFreq' epochs.
 */
static void TRAIN(times)
 int times;
{
 int i, report;

 report=ReportFreq;

 for (i=0; i<times; i++)
 {
  if (Epoch % report==0)     /* Time to report status */
  {

   /* If a single output unit, calculate the c index, if specified -dbr */
#if CALC_CINDEX
   if (Noutputs==1)
   {
    cTrain=C_INDEX(FirstOutput, NTrainingPatterns, TrainingInputs,TrainingOutputs);

    TrainSSE=TotalError;
    TrainErrorBits=TotalErrorBits;

    if (TestReporting==1)
    {
     cTest=C_INDEX(FirstOutput, NTestingPatterns, TestInputs,TestOutputs);

     TestSSE=TotalError;
     TestErrorBits=TotalErrorBits;
    }
   }
#endif /* CALC_CINDEX */

   if (!CALC_CINDEX || Noutputs>1)
   {
    TRAIN_SUMMARY();

    /* save "best" runs according c index or to RMS error  -phg */
    if (TestReporting==1)
    {
     TEST_SUMMARY();
    }
   }

   if (TestReporting==1) if (i>0 || (Epoch==EpochSaved && Epoch>0))
    RESHUFFLE_10BEST();

   if (Epoch==0 || Epoch==EpochSaved || i>0)
   {
    /* report out the results of this epoch's training (and testing) */
#define printargs\
"*--------------------------------------------------------------------*\n"
    printboth; 
#undef printargs
    if (Epoch==0 || Epoch==EpochSaved)
    {
#define printargs\
     "Epoch %d: \n", Epoch
     printboth; 
#undef printargs
    }
    else
    {
#define printargs\
     "Epoch %d: Did QuickProp on %6.2f %%, GradDesc on %6.2f %% of weights.\n",\
     Epoch, 100.0*(1.0 - (float)DidGradient/WeightCount), 100.0*(float)DidGradient/WeightCount
     printboth; 
#undef printargs
    }
    
    /* display the results */
    if (NTrainingPatterns>0)
    PRINT_RESULTS("TRAINING", NTrainingPatterns, TrainErrorBits, TrainSSE, (float)cTrain);

    if (TestReporting==1 && NTestingPatterns>0)
     PRINT_RESULTS(" TESTING", NTestingPatterns, TestErrorBits, TestSSE, (float)cTest);
   }

   if (!STOPPED_TRAINING() && i<times-1) 
   {
    KeepScore=1;

    TRAIN_ONE_EPOCH();

    KeepScore=0;  /* Turn off unnec. calcs now */
   }
   else i=times;
  }
  else  if (i<times-1)
  {
   TRAIN_ONE_EPOCH();
  }
 }
}
/****** END OF TRAIN() ******/


#if CALC_CINDEX
/*
 *  C_INDEX():
 *  Calculate c index for a particular unit k.
 *  Note: no provision for freeing/shrinking storage. -dbr
 *  Side effect is to write TotalError & TotalErrorBits
 */
static float C_INDEX(k, nPatterns, inputs, outputs)
 int k; /* unit of interest */
 int nPatterns;
 float **inputs;
 float **outputs;
{
 static int nPatternsMax=-1;
 static float *response=NULL;  /* float response[nPatternsMax] */
 static float *truth=NULL;  /* float truth[nPatternsMax] */

 if (nPatterns > nPatternsMax)
 {
  /* Expand or create response array */
  nPatternsMax=nPatterns;
  destroy(&response); /* no effect if response==NULL */
  ncreate(&response, nPatternsMax);
  destroy(&truth); /* no effect if truth==NULL */
  ncreate(&truth, nPatternsMax);
 }

 /* Copy to truth array */
 { register int iPatt;
   for(iPatt=0; iPatt<nPatterns; iPatt++)
     truth[iPatt]=outputs[iPatt][k-FirstOutput];
 }   /*  truth[iPatt][k-FirstOutput] */

 PREDICT_1(k, nPatterns, inputs, response);  /* Calculate response array */

 /* now generate %correct and RMSErr info for this epoch */
 TotalError=0.;
 TotalErrorBits=0;
 {int iPatt;
  for(iPatt=0; iPatt<nPatterns; iPatt++)
  {
   register float dif = truth[iPatt] - response[iPatt];
   TotalError +=dif*dif;
   if (((float) fabs((double) dif)) >=ScoreThreshold)
    TotalErrorBits++;
  }
 }
 return cindF((long)nPatterns, response, truth);
}
/****** END OF C_INDEX() ******/


/* 
 *  PREDICT_1():
 *  Calculate output response of 1 unit (unit k) for a set of input patterns.
 *  response[iPattern] gets output response of unit k to input pattern
 *  inputs[iPattern][] -dbr.
 *  
 */
static void PREDICT_1(k, nPatterns, inputs, response)
 int k; /* unit of interest */
 int nPatterns;
 float **inputs;
 float *response;
{
  register int iPatt;

 for(iPatt=0; iPatt<nPatterns; iPatt++)
 {
  FORWARD_PASS(inputs[iPatt]);
  response[iPatt]=Outputs[k];
 }
}
/****** END OF PREDICT_1() ******/
#endif /* CALC_CINDEX */


/*
 *  TRAIN_SUMMARY():
 *  Do a forward pass on TEST data to generate error stats  -phg
 */
static void TRAIN_SUMMARY()
{
 int j, k;

 KeepScore=1;
 TotalError=0.0;
 TotalErrorBits=0;

 for(j=0; j<NTestingPatterns; j++)
 {
  FORWARD_PASS(TrainingInputs[j]);

  for (k=FirstOutput; k<Nunits; k++)
  {
   (void) ERRFUN(TrainingOutputs[j][k-FirstOutput], Outputs[k]);
  }
 }
 TrainErrorBits=TotalErrorBits; /* Rename these for use later  -phg */
 TrainSSE=TotalError;
 KeepScore=0;
}


/*
 *  TEST_SUMMARY():
 *  Do a forward pass on TEST data to generate error stats  -phg
 */
static void TEST_SUMMARY()
{
 int j, k;

 KeepScore=1;
 TotalError=0.0;
 TotalErrorBits=0;

 for(j=0; j<NTestingPatterns; j++)
 {
  FORWARD_PASS(TestInputs[j]);

  for (k=FirstOutput; k<Nunits; k++)
  {
   (void) ERRFUN(TestOutputs[j][k-FirstOutput], Outputs[k]);
  }
 }
 TestErrorBits=TotalErrorBits; /* Rename these for use later  -phg */
 TestSSE=TotalError;
 KeepScore=0;
}
/****** END OF TEST_SUMMARY() ******/


/*
 *  RESHUFFLE_10BEST():
 *  Time to see if current prediction should be among cumulative best 10.
 */
static void RESHUFFLE_10BEST()
{
 int j, k, m, n;
 float cTesttemp=-1., RMSTesttemp=-1., CompValue;

 if (cTrain>=0.) cTesttemp=cTest;

 TestRMS=(float) sqrt((double)(TestSSE/(NTestingPatterns*Noutputs)));
 RMSTesttemp=(1.0 - TestRMS);

 
 for (j=0; j<10; j++) /* Compare current score to each of existing 10 */
 {
  if (BestByCindex && cTrain>=0.) CompValue=cTesttemp - BestcTest[j];
  else CompValue=RMSTesttemp - BestRMSTest[j];

  if (CompValue > 0) /* if current value is better than prev best -phg */
  {

   if (j==0) /* only if the very best is replaced, save the wts into an array */
   {
    for (m=0; m<Nunits; m++)
    {
     for (n=0; n<Nunits; n++)  BestWeights[m][n]=Weights[m][n];
    }
   }

   for (k=9; k>j; k--) /* Move others down one index from the one replaced */
   {
    BestEpoch[k]=BestEpoch[k-1];
    BestRMSTest[k]=BestRMSTest[k-1];
    if (cTrain>=0.) BestcTest[k]=BestcTest[k-1];
    ReportInterv[k]=ReportInterv[k-1];
   }

   /* Lastly, assign current index to new value  -phg */
   BestEpoch[j]=Epoch;  
   BestRMSTest[j]=(1.0 - TestRMS);
   if (cTrain>=0.) BestcTest[j]=cTest;
   ReportInterv[j]=ReportFreq;

   if (cTrain>=0.) cTesttemp=-1.0; /* reset */
   RMSTesttemp=-1.0;
  }
 }
}


/*
 *  TRAIN_ONE_EPOCH():
 *  Perform forward and back propagation once for each pattern in the
 *  training set, collecting deltas.  Then burn in the weights.
 *
 *  EpochWiseUpdate(per-epoch) vs per-pattern option offered.
 *  If using QuickProp, only EpochWiseUpdate available at present -
 *  check with Scott Fahlman for progress on this issue. -phg
 */
static void TRAIN_ONE_EPOCH()
{
 int i;

 DidGradient=0;  /* Reinitialize prior to training another epoch  -phg*/
 TotalError=0.0;
 TotalErrorBits=0;

 if (EpochWiseUpdate) CLEAR_SLOPES();

 for (i=0; i<NTrainingPatterns; i++)
 {
  if (!EpochWiseUpdate) CLEAR_SLOPES();

  FORWARD_PASS (TrainingInputs[i]);
  BACKWARD_PASS (TrainingOutputs[i]);

  if (!EpochWiseUpdate) UPDATE_WEIGHTS();
 }

 if (EpochWiseUpdate) UPDATE_WEIGHTS();

 Epoch++;
}
/****** END OF TRAIN_ONE_EPOCH() ******/


/*
 *  FORWARD_PASS():
 *  This is it, ya Habaayib:  the forward pass in BP.
 */
static void FORWARD_PASS(input)
 float input[];
{
 int i, j, *connect;
 float sum, *weight;

 /* Load in the input vector */
 for (i=0; i<Ninputs; i++)  Outputs[i+1]=input[i];

 /* For each unit, collect incoming activation and pass through sigmoid. */
 for (j=FirstHidden; j<Nunits; j++)
 {
  connect = Connections[j];
  weight = Weights[j];
  sum=0.0;
  i = Nconnections[j];
  while (i > 0)
  {
   i--;
   sum +=(Outputs[ connect[i] ] * weight[i]);
  }
  Outputs[j]=ACTIVATION(sum);
 }
}
/****** END OF FORWARD_PASS() ******/


/*
 *  ACTIVATION():
 *  Given the sum of weighted inputs, compute the unit's activation value.
 *  Defined unit types are SIGMOID and ASYMSIGMOID.
 */
static float ACTIVATION(sum)
 float sum;
{
 switch(Unit_type)
 {
  case SIGMOID:
   if (sum < -30.0)  return(-0.5);  /* -30 used to prevent overflow */
   return(float)(1.0/(1.0+(exp((double) -sum))-0.5));

  default:
   dbassert(Unit_type==ASYMSIGMOID);
   /* asymmetrical sigmoid function in range 0.0 to 1.0. */
   if (sum < -30.0)  return(0.0);  /* -30 used to prevent overflow */
   return(float)(1.0/(1.0+exp((double) -sum)));
 }
}
/****** END OF ACTIVATION() ******/


/*
 *  ERRFUN():
 *  Compute the error for one output unit.
 *  For reporting periods, calculate SSE and error bits.
 *  HyperErr==0 => use squared error.
 *  HyperErr==1 => use tanh.
 */
static float ERRFUN(desired, actual)
 float desired, actual;
{
 float dif;

 dif=desired - actual;

 if (KeepScore)
 {
  TotalError +=dif*dif;
  if (((float) fabs((double) dif)) >=ScoreThreshold)
  TotalErrorBits++;
 }

 if (HyperErr==0)  return (dif); /* Not using atanh for error */
 else  return ((float) log ((double) (1.0+dif) / (1.0-dif)));
}
/****** END OF ERRFUN() ******/


/*
 *  STOPPED_TRAINING():
 *  Stop training if beyond specified Epoch past peak generalization.
 */
static int STOPPED_TRAINING()
{
 if (SaveBestWts == 1 && Epoch>(MinEpochs + EpochSaved) && Epoch>(BeyondBestEpoch*BestEpoch[0]) && ForceTraining==0)
 {
  fflush(stdout);

#define printargs\
  "\nOOO-AHHH! STOPPED TRAINING: Epoch %d is at least %G X BestEpoch %d \n",\
  Epoch, BeyondBestEpoch, BestEpoch[0]
  printboth;
#undef printargs

  return(1);
 }
 else  return(0);
}


/*
 *  CLEAR_SLOPES():
 *  Save the current slope values as PrevSlopes, and "clear" all current
 *  slopes (actually set to corresponding weight, decayed a bit).
 */
static void CLEAR_SLOPES()
{
 int i, j;

 for (i=FirstHidden; i<Nunits; i++)
 {
  for (j=0; j<Nconnections[i]; j++)
  {
   PrevSlopes[i][j]=Slopes[i][j];
   Slopes[i][j]=(Decay * Weights[i][j]);
  }
 }
}
/****** END OF CLEAR_SLOPES() ******/


/*
 *  BACKWARD_PASS():
 *  Goal is a vector of desired values for the output units.  Propagate
 *  the error back through the net, accumulating weight deltas.
 */
static void BACKWARD_PASS(goal)
 float goal[];
{
 int i, j, cix, *connect;     /* cix is "connection index" */
 float errors, *slope, *weight;

 /* Compute error sums for output and other units  */
 for (i=FirstHidden; i<FirstOutput; i++)
  ErrorSums[i]=0.0;
 for (; i<Nunits; i++)
  ErrorSums[i]=ERRFUN(goal[i-FirstOutput], Outputs[i]);

 /* Back-prop.  When we reach a given unit in loop, error from all later
  * units will have been collected. */
 for (j=Nunits-1; j>=FirstHidden; j--)
 {
  connect = Connections[j];
  weight = Weights[j];
  slope = Slopes[j];
  errors=ACTIVATION_PRIME(Outputs[j]) * ErrorSums[j];
  i = Nconnections[j];
  while (i > 0)
  {
   i--;
   cix=connect[i];
   ErrorSums[cix] +=(errors * weight[i]);
   slope[i] +=(errors * Outputs[cix]);
  }
 }
}
/****** END OF BACKWARD_PASS() ******/


/*
 *  ACTIVATION_PRIME():
 *  Given the unit's activation value and sum of weighted inputs,
 *  compute the derivative of the activation with respect to the sum.
 *  Defined types are SIGMOID (-1 to +1) and ASYMSIGMOID (0 to +1).
 */
static float ACTIVATION_PRIME(value) /* modified for default return  -dbr */
 float value;
{
 switch(Unit_type)
 {
  case SIGMOID:
  /* Symmetrical sigmoid function. */
   return (SigmoidPrimeOffset + (0.25 -  value*value));

  default:
   dbassert(Unit_type==ASYMSIGMOID);
   /* asymmetrical sigmoid function in range 0.0 to 1.0. */
   return (SigmoidPrimeOffset + (value * (1.0 - value)));
 }
}
/****** END OF ACTIVATION_PRIME() ******/


/*
 *  UPDATE_WEIGHTS():
 *  Update all weights in the network as a function of each weight's current
 *  slope, previous slope, and the size of the last jump.
 */
static void UPDATE_WEIGHTS()
{
 int i, j;
 float next_step, shrink_factor, epsilon;

 shrink_factor=MaxFactor / (1.0 + MaxFactor);

 WeightCount=0;
 DidGradient=0;

 for (j=FirstHidden; j<Nunits; j++)
 {
  if (SplitEpsilon)
   epsilon = Epsilon / (float)Nconnections[j];
  else
   epsilon = Epsilon;
  for (i=0; i<Nconnections[j]; i++)
  {
   next_step=0.0;
   WeightCount++;
   if ((DeltaWeights[j][i] > ModeSwitchThreshold) && UseQuickProp==1)
   {                            /* Last step was signif. +ive..... */
    if (Slopes[j][i] > 0.0)  /* Add in epsilon if +ive slope */
     next_step +=(epsilon * Slopes[j][i]);
    /* If slope > (or close to) prev slope, take max size step. */
    if (Slopes[j][i] > (shrink_factor * PrevSlopes[j][i]))
     next_step +=(MaxFactor * DeltaWeights[j][i]);
    else        /*  Use quadratic estimate */
     next_step +=((Slopes[j][i]/(PrevSlopes[j][i]-Slopes[j][i]))
                     * DeltaWeights[j][i]);
   }
   else if ((DeltaWeights[j][i] < -ModeSwitchThreshold) && UseQuickProp==1)
   {                          /* Last step was signif. -ive.... */
    if (Slopes[j][i] < 0.0)/* Add in epsilon if -ive slope */
     next_step +=(epsilon * Slopes[j][i]);
    /* If slope < (or close to) prev slope, take max size step. */
    if (Slopes[j][i] < (shrink_factor * PrevSlopes[j][i]))
     next_step +=(MaxFactor * DeltaWeights[j][i]);
    else        /*  Use quadratic estimate */
     next_step +=((Slopes[j][i]/(PrevSlopes[j][i]-Slopes[j][i]))
                     * DeltaWeights[j][i]);
   }
   else       /* Normal gradient descent, complete with momentum */
   {
    DidGradient++;
    next_step +=((epsilon * Slopes[j][i])
       + (Momentum * DeltaWeights[j][i]));
   }

   /* Set delta weight, and adjust the weight itself. */
   DeltaWeights[j][i]=next_step;
   Weights[j][i] +=next_step;
  }
 }
}
/****** END OF UPDATE_WEIGHTS() ******/


/*
 *  PRINT_RESULTS():
 *  Display line with error information.
 */
static void PRINT_RESULTS(label, nPatterns, errorBits, SSE, cIndex)
 char *label;
 int nPatterns;
 int errorBits;
 float SSE;
 float cIndex; /* Negative if unknown or undefined */
{
#define printargs\
 " %s: %6.2f %%correct ;  RMSErr=%7.5f", label,\
          (float) 100*(1-(float)errorBits/(nPatterns*Noutputs)),\
          (float) sqrt((double)(SSE/(nPatterns*Noutputs)))
 printboth; 
#undef printargs
 if (cIndex>=0.)
#define printargs\
 " ;  c index=%7.5f", cIndex
 printboth; 
#undef printargs
#define printargs\
 "\n"
 printboth; 
#undef printargs
}
/****** END OF PRINT_RESULTS() ******/


/*
 *  PRINT_10BEST():
 *  Display the best epochs of training.
 */
static void PRINT_10BEST()
{
 int k, m, n, q;
 int BestEpochtemp;
 float RelScore[10], MinTest, Crange, RMSrange;

 /* If reporting on test dataset, sort 10 best generalizing epochs */
 /* This sorting is needed to display graphic by epoch no. */
 if (TestReporting==1)
 {
  for (n=0; n<10; n++)  /* Run each Best value by */
  {
   BestEpochtemp=BestEpoch[n];
   for (k=0; k<10; k++)  /* Sort each by decreasing Epoch */
   {
    if (BestEpochtemp > SortedEpoch[k])
    {
     for (m=9; m>k; m--) /* Move each down to proper sorted level */
     {
      SortedEpoch[m]=SortedEpoch[m-1];
      SortedRMSTest[m]=SortedRMSTest[m-1];
      if (cTrain>=0.) SortedcTest[m]=SortedcTest[m-1];
      SortedReportInterv[m]=SortedReportInterv[m-1];
     }
     SortedEpoch[k]=BestEpoch[n];
     SortedRMSTest[k]=BestRMSTest[n];
     if (cTrain>=0.) SortedcTest[k]=BestcTest[n];
     SortedReportInterv[k]=ReportInterv[n];
     BestEpochtemp=-1;
    }
   }
  }

  /* Generate relative score to create graphic "***" output */
  if (BestByCindex && cTrain>=0.)
  {
   /* If <10 scores recorded, use only real value as min */
   MinTest=BestcTest[9];
   for (q=9; q>0; q--)
   {
    if (BestcTest[q] < 0. && BestEpoch[q-1]!=0) MinTest=(BestcTest[q-1]);
    if (BestcTest[q] < 0. && BestEpoch[q-1]==0) MinTest=(BestcTest[q-2]);
   }

   /* if there's no significant variation, set all to 10 asterixes */
   if ((Crange=(BestcTest[0]-MinTest)) < 0.000001)
    for (q=0; q<10; q++) RelScore[q]=10;

   /* otherwise, give 1 to 10 asterixes according to relative range */
   else for (q=0; q<10; q++)
   {
    RelScore[q]=(SortedcTest[q]-MinTest)/Crange;
    RelScore[q]=RelScore[q]*9 + 1; /* assures min 1 and max of 10 */
   }
  }
  else
  {
   /* If <10 scores recorded, use only real value as min */
   MinTest=BestRMSTest[9];
   for (q=9; q>0; q--)
   {
    if (BestRMSTest[q] < 0. && BestEpoch[q-1]!=0) MinTest=(BestRMSTest[q-1]);
    if (BestRMSTest[q] < 0. && BestEpoch[q-1]==0) MinTest=(BestRMSTest[q-2]);
   }
   if ((RMSrange=(BestRMSTest[0]-MinTest)) < 0.000001)
    for (q=0; q<10; q++) RelScore[q]=10;
   else for (q=0; q<10; q++)
   {
    RelScore[q]=(SortedRMSTest[q]-MinTest)/RMSrange;
    RelScore[q]=RelScore[q]*9 + 1; /* assures min 1 and max of 10 */
   }
  }
  if (SaveBestWts == 1)  PRINT_GRAPHICS(stdout, RelScore);
  if (resfile && SaveBestWts == 1) PRINT_GRAPHICS(resfile, RelScore);
 }
}
/****** END OF PRINT_10BEST() ******/


/*
 *  PRINT_GRAPHICS():
 *  Display horizontal bar graph of up to ten best generalizing reported epochs.
 */
static void PRINT_GRAPHICS(stream, RelScore)
 FILE *stream;
 float RelScore[10];
{
 int i, o;
 int SkipSome, EpochPast;

 fprintf(stream,"\n\
*--------------------------------------------------------------------*\n");

#if CALC_CINDEX
  fprintf(stream,"\n(Up to) 10 BEST runs on the TEST set, according to");
  if (BestByCindex) fprintf(stream,"\
 c index:\n");
  else              fprintf(stream,"\
 RMS error:\n");
 if (ReportRandSeed==1)  fprintf(stream,"\
  (using training weight SEED=%lu)\n", seed);

  fprintf(stream,"                                                  ");
  fprintf(stream,"  RELATIVE\n");

  if (BestByCindex)
  {
   fprintf(stream,"                                                  ");
   fprintf(stream,"  C INDEX\n");
  }
  else
  {
   fprintf(stream,"                                                  ");
   fprintf(stream," (1-RMSErr)\n");
  }
   fprintf(stream,"                                                  ");
   fprintf(stream," LOW-----HI\n\n");
#endif /* CALC_CINDEX */

#if !CALC_CINDEX
   fprintf(stream,"\n10 BEST runs on the TEST set, according to RMS error:\n");
   fprintf(stream,"  (using training weight SEED=%lu)\n", seed);
   fprintf(stream,"\n                              ");
   fprintf(stream,"  RELATIVE\n");
   fprintf(stream,"                              ");
   fprintf(stream," (1-RMSErr)\n");
   fprintf(stream,"                              ");
   fprintf(stream," LOW-----HI\n\n");
#endif /* CALC_CINDEX */

   EpochPast=0;
   SkipSome=0;

   for (i=9; i>=0; i--)
   {
    if (SortedEpoch[i] > 0)
    {
     if ((i!=9) && (SortedEpoch[i]!=EpochPast+SortedReportInterv[i]) && SortedEpoch[i]>EpochSaved)
     {
      fprintf(stream,"\
@.................................................................\n");
      SkipSome++;
     }
     fprintf(stream," Epoch:%6d   RMSErr=%8.6f ", SortedEpoch[i], (1.0 - SortedRMSTest[i]));
     if (cTrain>=0.) fprintf(stream,"  c index=%8.6f ", SortedcTest[i]);

     for (o=0; o<(int)RelScore[i]; o++) fprintf(stream,"*");
     fprintf(stream,"\n");
    }
    if (SortedEpoch[i] >= 0) EpochPast=SortedEpoch[i];
   }
   if (SkipSome > 0) fprintf(stream,"\n\
@..............................=> Skips at least one reporting period.");
   fprintf(stream,"\nTraining halted after Epoch %d.", Epoch);
   fprintf(stream,"\n*--------------------------------------------------------------------*\n");
}
/****** END OF PRINT_GRAPHICS() ******/


/*
 *  TEST_SAVE():
 *  Save predictions to a .prd file -- depending on settings, either w/o training,
 *  or based on the final trained weights orbest generalizing weights found.
 */
static void TEST_SAVE()
{
 FILE  *outfile, *fopen();
 int  i, j, m, n, response=0;
 char realfname[256], fprdname[256];

 if (IBM==0)   /* Use long filename format if not an IBM or DOS machine  -phg */
  sprintf(realfname, "%s-%d-%d-%d.prd", fname, Ninputs, Nhidden, Noutputs);
 else
  sprintf(realfname, "%s.prd", fname);

 outfile=fopen (realfname, "r");
 if (!outfile) outfile=fopen (realfname, "w"); 
 else 
 {
  fclose (outfile);
  if (PredictOnly==1) fprintf(stderr,"\
\n*--------------------------------------------------------------------*\n");
  fprintf(stderr," >Prediction file \"%s\" already exists, OVERWRITE?[y] ", realfname);
  if (response=READ_DEFAULT(1)) outfile=fopen (realfname, "w");
  else
  {
   fprintf(stderr,"  >Enter new name for prediction file [assumes .prd]: ");
   gets(fprdname);
   sprintf(realfname, "%s.prd", fprdname);
   outfile=fopen (realfname, "w");
  }
 }

 fprintf(outfile, "SEQU"); /* Sequence number as in ID in header row  -phg */
 for (i=0; i<Noutputs; i++)
 {
  fprintf(outfile," PRED%d",i+1); /* Then, PRED#s in header row */
 }
 for (i=0; i<Noutputs; i++)
 {
  fprintf(outfile," TRUE%d",i+1); /* Then, TRUE#s in header row */
 }
  fprintf(outfile,"\n");

 if (SaveBestWts)
 {
  for (m=0; m<Nunits; m++)
   for (n=0; n<Nunits; n++)  Weights[m][n]=BestWeights[m][n];
 }

 for(j=0; j<NTestingPatterns; j++)
 {
  FORWARD_PASS(TestInputs[j]);
  fprintf(outfile, "%d",j+1);
  for (i=0; i<Noutputs; i++)
   fprintf(outfile, " %f", Outputs[FirstOutput+i]);
  for (i=FirstOutput; i<Nunits; i++)
   fprintf(outfile, " %f", TestOutputs[j][i-FirstOutput]);
   fprintf(outfile, "\n");
 }
 fclose (outfile);
 
 if (response)
 {
#define printargs\
  "  ... Predictions overwritten on existing \"%s\"\n\n", realfname
  printboth;
#undef printargs
 }
 else
 {
#define printargs\
  "   ... Predictions saved to \"%s\"\n\n", realfname
  printboth;
#undef printargs
 }
}
/****** END OF TEST_SAVE() ******/


/*
 *  DUMP_WEIGHTS():
 *  Save weights to a .wts file -- depending on settings, either
 *  the final trained weights, or best generalizing weights found. -phg
 */
static void DUMP_WEIGHTS()
{
 FILE  *outfile, *fopen();
 int  i, j, m, n, response=0;
 char realfname[256], fwtname[256];

 if (IBM==0)   /* Use long filename format if not an IBM or DOS machine  -phg */
  sprintf(realfname, "%s-%d-%d-%d.wts", fname, Ninputs, Nhidden, Noutputs);
 else
  sprintf(realfname, "%s.wts", fname);

 outfile=fopen (realfname, "r");
 if (!outfile) outfile=fopen (realfname, "w");
 else 
 {
  fclose (outfile);
  fflush(stdout);
  fprintf(stderr," >Weights file \"%s\" already exists, OVERWRITE?[y] ", realfname);
  if (response=READ_DEFAULT(1)) outfile=fopen (realfname, "w");
  else
  {
   fflush(stdout);
   fprintf(stderr,"  >Enter new name for weights file [assumes .wts]: ");
   gets(fwtname);
   sprintf(realfname, "%s.wts", fwtname);
   outfile=fopen (realfname, "w");
  }
 }
 
 EpochSaved = Epoch;

 if (SaveBestWts)
 {
  EpochSaved = BestEpoch[0];
  for (m=0; m<Nunits; m++)
   for (n=0; n<Nunits; n++)  Weights[m][n]=BestWeights[m][n];
 }

 /* Indicate Epoch at which wts were frozen */
 fprintf(outfile, "EpochSaved %d\n", EpochSaved);

 for (j=FirstHidden; j<Nunits; j++)
  for (i=0; i<Nconnections[j]; i++)
   /* Put each connection and weight on unique line  -phg */
   fprintf(outfile, "%d %d %g\n", j, Connections[j][i], Weights[j][i]);

 fclose (outfile);
   
 if (response)
 {
#define printargs\
  "  ... Weights overwritten on existing \"%s\"\n", realfname
  printboth;
#undef printargs
 }
 else
 {
#define printargs\
  "   ... Weights saved to \"%s\"\n", realfname
  printboth;
#undef printargs
 }
}
/****** END OF DUMP_WEIGHTS() ******/


/*
 *  FREE_MEMORY():
 *  Do it.
 */
static void FREE_MEMORY()
{
 int i;

 for(i=0;i<NTrainingPatterns;i++)
 {
  free((void *)TrainingInputs[i]);
  free((void *)TrainingOutputs[i]);
 }
 free((void *)TrainingInputs);
 free((void *)TrainingOutputs);

 for(i=0;i<NTestingPatterns;i++)
 {
  free((void *)TestInputs[i]);
  free((void *)TestOutputs[i]);
 }
 free((void *)TestInputs);
 free((void *)TestOutputs);

 for(i=0;i<Nunits;i++)
 {
  free((void *)Connections[i]);
  free((void *)Weights[i]);
  free((void *)BestWeights[i]);
  free((void *)DeltaWeights[i]);
  free((void *)Slopes[i]);
  free((void *)PrevSlopes[i]);
 }
 free((void *)Connections);
 free((void *)Weights);
 free((void *)BestWeights);
 free((void *)Slopes);
 free((void *)PrevSlopes);
 free((void *)DeltaWeights);

 free((void *)Outputs);
 free((void *)ErrorSums);
 free((void *)Nconnections);
}


/* * * * * * * * * * * * * * END of NevProp C code * * * * * * * * * * */
@EOF

chmod 640 NevProp/np.c

echo x - NevProp/cind.c
cat >NevProp/cind.c <<'@EOF'
/* Author: David B. Rosen
   Copyright (C) 1993 by the author and/or University of Nevada.
   All rights reserved. 

Change log:

10 May 1993  DBR  Remember to destroy() before exceptional return
04 Mar 1993  DBR  Return 3 distinct codes defined in .h.  Don't crash.

*/

/* system: */
#include "stdio.h"  /* SunOS 4.1.3 assert.h uses stderr but doesn't define! */
#include "assert.h"

/* general utility: */
#include "const.h" /* define const if necessary */
#include "do.h"
#include "create.h" /* needs create.c */

/* other: */
#include "rank.h" /* needs rank.c */

#include "cind.h" /* this */

float cindF(n, predict, actual)
     long n;
     const float *predict;
     const float *actual;
{
  float *rank;
  register long n1=0; /* will be number of actual==hi points */
  register float sum=0.0;  /* will be sum of ranks of actual==hi points */
  register float hi,lo;  /* the two possible values of actual */

  if (n==0) return CIND_ZERO_ACTUAL_VALS;

  hi=actual[0]; lo=actual[0];
  do0L(i, n,
       if (actual[i] != hi)
       { if (actual[i] < hi)
	   lo = actual[i];
         else
	   hi=actual[i];
	 break;
       }
  );

  if (hi==lo) return CIND_ONE_ACTUAL_VAL;

  /* create and calculate rank array: */
  ncreate(&rank, n);
  rankF(n, predict, rank);

  /* calculate sum, n1 */
  do0(i, n,
    if (actual[i] == hi)
    { n1++;
      sum += rank[i];
    } else if (actual[i]!=lo)
    { destroy(&rank);
      return CIND_OVER_TWO_ACTUAL_VALS;
    }
  );

  destroy(&rank);

  assert((n1 != 0) && (n-n1 != 0)); /* shouldn't happen because we already
				    found distinct hi and lo above */

  { register float mean = sum/n1;  /* mean rank of actual==hi points */
    return (mean - (n1 + 1)*.5 ) / (n - n1) ;
  }
}
@EOF

chmod 640 NevProp/cind.c

echo x - NevProp/create.c
cat >NevProp/create.c <<'@EOF'
/* Copyright (C) 1991,1993 by David B. Rosen.  All rights reserved. */

#include "stdio.h"
#include "stdlib.h"

#include "create.h"

void _ncreate(ptrptr,size)
     void **ptrptr;
     unsigned long size; /* avoid size_t for pre-ansi  -DBR 5/93 */
{
    if (!(*ptrptr = (void *)malloc(size))) 
      {
	  fprintf(stderr,
	  "Memory allocation failure for %d bytes in _ncreate().  Exiting.\n",
		  (int)size);
	  exit(1);
      }
}

void _destroy(ptrptr)
     void **ptrptr;
{
  if (*ptrptr) free(*ptrptr);
  *ptrptr = (void *)0;
}
@EOF

chmod 640 NevProp/create.c

echo x - NevProp/rank.c
cat >NevProp/rank.c <<'@EOF'
/* Author: David B. Rosen
   Copyright (C) 1993 by the author and/or University of Nevada.
   All rights reserved. */

#include "stdlib.h"

#include "const.h" /* define const if necessary */
#include "do.h"
#include "create.h"

#include "rank.h"

/* value-index type: */
typedef struct {float value; long index;} viT;

static int compare(a, b)
     const viT *a;
     const viT *b;
{
  if (a->value < b->value)
     return -1;
  else if (a->value > b->value)
    return 1;
  else
    return 0;
}

void rankF(n, array, rank)
     long n;
     const float *array;  /*  array[n]  */
     float *rank;  /* rank[n]  */
{
  viT *vi;

  /* create value-index array: */
  ncreate(&vi, n);
  do0L(i, n,
      vi[i].value = array[i];
      vi[i].index = i;
  );

  /* sort the v i array: */
  qsort(vi, n, sizeof(vi[0]), 
	(int (*)proto((const void *, const void *)))compare);

  /* write the rank: */
  { register long iFirst=0;
    while (iFirst < n) /* loop over groups of consecutive tied values */
    { register long iLast=iFirst;
      while(iLast+1 < n) /* loop over tied values in this group */
      { if (vi[iLast+1].value != vi[iFirst].value) break;
	iLast++;
      }
      { register float thisRank=.5*(iLast+iFirst+2);
	doup(i, iFirst, iLast, rank[vi[i].index] = thisRank;  );
      }
      iFirst=iLast+1;
    }

    destroy(&vi);
  }
}
@EOF

chmod 640 NevProp/rank.c

echo x - NevProp/cind.h
cat >NevProp/cind.h <<'@EOF'
#ifndef _CIND_H
#define _CIND_H

/* cindF -- C INDex using Float.
   Returns C Index of (predict[i],actual[i]), i = 0..n-1.  No side effects. 
   Each value actual[i] should have one of *two* distinct values (one "high",
   the other "low").   Returns a negative (or perhaps NaN-valued on some
   systems) constant if c index is undefined for various reasons.  (See names
   of constants below.)

*/

#include "const.h" /* define const if necessary */
#include "protoTyp.h"

extern float cindF proto((long n, const float *predict, const float *actual));

#define CIND_OVER_TWO_ACTUAL_VALS (-991.0)  /* too many distinct values */
#define CIND_ONE_ACTUAL_VAL (-992.0)  /* only one distinct value */
#define CIND_ZERO_ACTUAL_VALS (-993.0)  /* no values (n == 0) */
/* numeric values above subject to change of course! */

/*
Change log:

04 Mar 1993  DBR  Return 3 distinct exception constants;  tolerate and return
                  CIND_OVER_TWO_ACTUAL_VALS.
*/

#endif /* _CIND_H */
@EOF

chmod 640 NevProp/cind.h

echo x - NevProp/const.h
cat >NevProp/const.h <<'@EOF'
/* define const as a null token for non-ansi C */
#ifndef __STDC__

#ifndef const
#define const
#endif

#endif
@EOF

chmod 640 NevProp/const.h

echo x - NevProp/create.h
cat >NevProp/create.h <<'@EOF'
#ifndef _CREATE_H
#define _CREATE_H
/* 

Copyright (C) 1991, 1992, 1993 by David B. Rosen.  All rights reserved.
Last Modified: 2 May 1993

type *ptr;

create(&ptr)
   Allocates a single object of the type of *ptr and sets ptr to point to it.

ncreate(&ptr,n) 
   Allocates an array of n objects of the type of *ptr and sets ptr to point 
   to 'em.

createvar(&ptr, fieldName, fieldElements)
   Allocates a single object of the type of *ptr and sets ptr to point to it.
   The object has a variable-length final field fieldName[], whose allocated
   length will be fieldElements.

destroy(&ptr)
   Deallocates memory allocated by the allocation functions above.
   and sets ptr = NULL.   ptr may already be NULL, in which case destroy()
   simply returns and has no effect.

*/

#include "protoTyp.h"

void _ncreate proto((void **, unsigned long)); /* avoid needing size_t */
void _destroy proto((void **));

#define create(ptrptr) ncreate((ptrptr),1)
#define ncreate(ptrptr,n) _ncreate((void **)(ptrptr), \
 (unsigned long)((n)*sizeof(**(ptrptr))))
#define createvar(ptrptr,fieldName,fieldElements) \
 _ncreate((void **)(ptrptr), (unsigned long)(sizeof(**(ptrptr)) + \
 (fieldElements)*sizeof((**(ptrptr)).fieldName[0]) ))

#define destroy(ptrptr) _destroy((void **)(ptrptr));

#endif /* _CREATE_H */
@EOF

chmod 640 NevProp/create.h

echo x - NevProp/do.h
cat >NevProp/do.h <<'@EOF'
#ifndef _DO_H
#define _DO_H

/* Copyright (C) 1991, 1992, 1993 by David B. Rosen.  All rights reserved. */

/*

do.h (See change log below for version number.)

  nrepeat(n, body)
   is a statement that executes body n times serially. Note that n
   is simply an expression that tells nrepeat how many times to 
   execute the body.  There is no loop variable available within the 
   loop. If you need to use a loop variable, use one of the macros 
   below instead. 

  doup(ivar,n1,n2,body)
   is a statement that executes body serially, with ivar
   successively taking on the ascending values n1, n1+1, ..., n2.

  dodown(ivar,n1,n2,body)
   is a statement that executes body serially, with ivar
   successively taking on the descending values n1, n1-1, ..., n2.

  do0(ivar,n,body)
   is a statement that executes body n times serially, with ivar
   successively taking on the descending values n-1, n-2, ... 0.   If
   you want the values of ivar to ascend, use doup instead.

  do1(ivar,n,body)
   is a statement that executes body n times serially, with ivar
   successively taking on the descending values n, n-1, n-2, ... 1.  If
   you want the values of ivar to ascend, use doup instead.

For all of the loop macros described above:
   o  The loop variable you name (ivar above) is defined by the macro as a
       register int local to the loop.  Do not declare it yourself. Do not
       modify its value, as the result of doing so should be considered
       unpredictable.
   o  The macro itself is syntactically a single statement. This means you
       should normally follow it by a semicolon, as you would any statement
       that does not end with ``}''.
   o  LIMITATION:  The body may not contain any commas except enclosed within
       parentheses () or within a string literal, otherwise a(n ANSI, but 
       perhaps not a traditional,) C preprocessor will think these delimit
       additional (spurious) arguments for the loop macro!
   o  The count (n) or range of values (n1 and n2) are expressions (rvalues)
       that are each evaluated only once,  before entering the loop the
       first time.  They can be of any type that can be implicitly converted
       to int, but do not use an expression that may be out of range of a
       (signed) int.
   o  break and continue may be used in the normal way within these loops.
   o  A final semicolon is automatically appended to the ``body'', so you need
       not supply one yourself, though you may if you wish.
   o  The body is enclosed directly in {} by the macro, so the body may, e.g.,
       begin with local declarations and definitions. Because you can't use
       commas (see above), use e.g. ``int i; int j;'' instead of ``int i, j;''.
   o  Problem:  The source line number reported in case of error
       anywhere within the body may reflect the first line of the macro
       rather than the line where the error occurred, depending on your
       preprocessor and compiler.
   o  Problem: C-specific editor modes and pretty-print software will 
       probably indent the loop as a function with arguments rather than 
       as a control structure with a body.

Note: do0 and do1 (and nrepeat) are likely to run slightly faster
than doup or dodown.  This is because doup and dodown have to
compare ivar to n2 every time through the loop, while do0 and do1 need
only test the sign of ivar.

For each of the macro names above, appending the letter `L' (do0L,
nrepeatL, etc.) gives a version of the macro that uses a long loop
variable instead of int.

*/

#if 0 /* example, so it can contain comments */

Example:
  #define N 100
  int a[N];
  ...
  ...
  nrepeat(5,
    /* body to be repeated 5 times */
    do0(i,N,    /* note that we needn't declare i */
      /* body to be repeated with i = N-1, N-2, ..., 2, 1, 0 */
      if (a[i] > 0 )
        func1(a[i]);
      else 
        func2(a[i]);
      doup(j, i, N-1, func3(a[i], a[j]) ); /* we chose to omit `;' in body */
    ); /* end of do0 i */
  ); /* end of nrepeat */

#endif /* example */

/*
To do:

  * implement do0L, nrepeatL, etc. using long loop variable.
  (currently only implemented on machines that happen to have int==long.)

*/

/* End of external specification and user documentation. 
   Implementation follows. */

/* 

Note: each macro is enclosed in a do {...} while(0) to make it syntactically
   a single statement, so as not to screw up if statements, etc. 

 Change log:

   Version 1.02:

   16 Mar 1993: Implemented L-suffix (long) portably.

   Version 1.01:

   02 Feb 1993: Added long (L-suffix) versions to documentation;
   implemented for int==long machines only.

   Version 1.0:

   18 Nov 1992 DBR - Improved documentation above.


   27 Oct 1992 DBR - Isolated and documented comma problem.  
                     Fixed buggy comma within doup and dodown macro defs.
   Version 0.8:

   1991
*/

#define do1(ivar,n,body) _do1Type(int,ivar,n,body)
#define do0(ivar,n,body) _do0Type(int,ivar,n,body)
#define doup(ivar,n1,n2,body) _doupType(int,ivar,n1,n2,body)
#define dodown(ivar,n1,n2,body) _dodownType(int,ivar,n1,n2,body)
#define nrepeat(n,body) _nrepeatType(int,n,body)

#define do1L(ivar,n,body) _do1Type(long,ivar,n,body)
#define do0L(ivar,n,body) _do0Type(long,ivar,n,body)
#define doupL(ivar,n1,n2,body) _doupType(long,ivar,n1,n2,body)
#define dodownL(ivar,n1,n2,body) _dodownType(long,ivar,n1,n2,body)
#define nrepeatL(n,body) _nrepeatType(long,n,body)

#define _do1Type(type,ivar,n,body) \
do { \
    register type (ivar)=((n)+1); \
    while ( --(ivar) > 0 ) \
      { \
	  body; /* include final ';' so user needn't for single-stmt body */ \
      } \
    } while(0)

#define _do0Type(type,ivar,n,body) \
do { \
    register type (ivar)=(n); \
    while ( --(ivar) >= 0 ) \
      { \
	  body; /* include final ';' so user needn't for single-stmt body */ \
      } \
    } while(0)

#define _doupType(type,ivar,n1,n2,body) \
do { \
    register type (ivar)=(n1)-1; register type _doup_final_=(n2); \
    while( (ivar)++ < _doup_final_ ) \
      { \
	  body; /* include final ';' so user needn't for single-stmt body */ \
      } \
    } while(0)

#define _dodownType(type,ivar,n1,n2,body) \
do { \
    register type (ivar)=(n1)+1; register type _doup_final_=(n2); \
    while( --(ivar) >= _doup_final_ ) \
      { \
	  body; /* include final ';' so user needn't for single-stmt body */ \
      } \
    } while(0)

#define _nrepeatType(type,n,body) _do0Type(type, _i_nrepeat_, (n), body)

#endif /* _DO_H */
@EOF

chmod 640 NevProp/do.h

echo x - NevProp/rank.h
cat >NevProp/rank.h <<'@EOF'
#ifndef _RANK_H
#define _RANK_H


#include "const.h" /* define const if necessary */
#include "protoTyp.h"

/* write rank[i] as rank (1-based) of array[i].  Ties give average rank
   among tied elements. 
   E.g. array = {3, 3, 9, 2.8} gives rank = {2.5, 2.5, 4, 1}*/ 
extern void rankF proto(( /* `F' for float */
		  long n
		  , const float *array  /*  array[n]  */
		  , float *rank  /* rank[n] (elements written) */
		  ));

#endif /* _RANK_H */
@EOF

chmod 640 NevProp/rank.h

echo x - NevProp/protoTyp.h
cat >NevProp/protoTyp.h <<'@EOF'
/* lets you declare functions as for example
       void myfunc proto((float **, int));
   and have this work in both ansi and traditional C

   Renamed from proto.h to protoArgs.h to avoid conflict with another
   proto.h distributed with gcc.
   Then renamed to protoTyp.h to be short enough for MS-DOS.

*/
   
#ifndef proto
#ifdef __STDC__
#define proto(a) a
#else
#define proto(a) ()
#endif /* __STDC__ */
#endif /* proto */
@EOF

chmod 640 NevProp/protoTyp.h

echo x - NevProp/iris.net
cat >NevProp/iris.net <<'@EOF'

******THIS IS A SAMPLE .net FILE-- see NP_READ.ME for more info*****

This is part of Fischer's iris database. The 50 members of species1 
have been eliminated, leaving 50 members each of species2 and species3. 
The flowers are listed in random order.

The four predictor variables are measurements of flower parts.
The outcome 'species2'=1 if the flower truly belongs to species2, 
and =0 if it belongs to species3. 

sl  sw   pl   lw     species

#

IBMorDOS 0            UseQuickProp 1        EpochWiseUpdate 1 

BestByCindex 0        MinEpochs   200       BeyondBestEpoch 2.0 

WeightRange 0.01      HyperErr 1            SigmoidPrimeOffset 0.1
Epsilon 0.1           SplitEpsilon 1        Momentum 0.1
WeightDecay -0.001    ScoreThreshold 0.1

QPMaxFactor 1.75      QPModeSwitchThreshold 0.0

Ninputs 4             Nhidden 3           Noutputs 1
UnitType 2
Connectcalls 2
1  4      5  7
5  7      8  8

NTrainingPatterns 50   NTestingPatterns 50

5.000 2.000 3.500 1.000    1
5.900 3.000 4.200 1.500    1
6.100 2.800 4.000 1.300    1
6.100 2.800 4.700 1.200    1
6.000 3.000 4.800 1.800    0
5.800 2.700 4.100 1.000    1
6.700 3.100 4.400 1.400    1
6.500 3.000 5.500 1.800    0
5.100 2.500 3.000 1.100    1
6.700 3.000 5.200 2.300    0
5.500 2.400 3.800 1.100    1
6.900 3.100 5.400 2.100    0
6.100 3.000 4.900 1.800    0
6.200 2.900 4.300 1.300    1
6.400 3.100 5.500 1.800    0
5.700 2.600 3.500 1.000    1
4.900 2.500 4.500 1.700    0
5.500 2.600 4.400 1.200    1
7.700 2.600 6.900 2.300    0
7.600 3.000 6.600 2.100    0
6.100 2.900 4.700 1.400    1
6.300 2.500 4.900 1.500    1
5.700 2.800 4.100 1.300    1
5.700 3.000 4.200 1.200    1
5.500 2.400 3.700 1.000    1
5.600 2.700 4.200 1.300    1
5.700 2.800 4.500 1.300    1
6.200 2.200 4.500 1.500    1
6.900 3.200 5.700 2.300    0
6.900 3.100 5.100 2.300    0
6.800 2.800 4.800 1.400    1
6.300 3.300 6.000 2.500    0
7.200 3.200 6.000 1.800    0
6.600 3.000 4.400 1.400    1
6.700 3.100 4.700 1.500    1
6.400 2.800 5.600 2.200    0
7.700 3.800 6.700 2.200    0
6.400 3.200 5.300 2.300    0
6.500 3.000 5.800 2.200    0
5.600 2.800 4.900 2.000    0
5.600 3.000 4.500 1.500    1
5.700 2.500 5.000 2.000    0
6.500 2.800 4.600 1.500    1
6.500 3.000 5.200 2.000    0
5.800 2.700 5.100 1.900    0
6.300 3.300 4.700 1.600    1
7.100 3.000 5.900 2.100    0
6.000 2.200 5.000 1.500    0
6.400 2.700 5.300 1.900    0
6.000 2.700 5.100 1.600    1

5.600 2.900 3.600 1.300    1
6.800 3.200 5.900 2.300    0
6.800 3.000 5.500 2.100    0
7.900 3.800 6.400 2.000    0
7.000 3.200 4.700 1.400    1
6.400 2.900 4.300 1.300    1
7.700 3.000 6.100 2.300    0
6.700 3.100 5.600 2.400    0
6.300 2.300 4.400 1.300    1
5.700 2.900 4.200 1.300    1
5.800 2.600 4.000 1.200    1
6.100 2.600 5.600 1.400    0
5.600 2.500 3.900 1.100    1
7.700 2.800 6.700 2.000    0
6.300 2.500 5.000 1.900    0
5.500 2.500 4.000 1.300    1
5.800 2.700 3.900 1.200    1
6.300 2.900 5.600 1.800    0
7.200 3.600 6.100 2.500    0
5.200 2.700 3.900 1.400    1
5.800 2.700 5.100 1.900    0
6.000 3.400 4.500 1.600    1
6.900 3.100 4.900 1.500    1
5.000 2.300 3.300 1.000    1
5.500 2.300 4.000 1.300    1
6.000 2.900 4.500 1.500    1
6.400 2.800 5.600 2.100    0
5.400 3.000 4.500 1.500    1
5.900 3.000 5.100 1.800    0
6.300 2.700 4.900 1.800    0
6.200 2.800 4.800 1.800    0
7.200 3.000 5.800 1.600    0
6.600 2.900 4.600 1.300    1
6.000 2.200 4.000 1.000    1
6.500 3.200 5.100 2.000    0
6.100 3.000 4.600 1.400    1
4.900 2.400 3.300 1.000    1
5.600 3.000 4.100 1.300    1
7.400 2.800 6.100 1.900    0
6.200 3.400 5.400 2.300    0
6.700 3.300 5.700 2.100    0
6.700 3.000 5.000 1.700    1
6.300 3.400 5.600 2.400    0
6.700 2.500 5.800 1.800    0
5.800 2.800 5.100 2.400    0
7.300 2.900 6.300 1.800    0
6.700 3.300 5.700 2.500    0
5.900 3.200 4.800 1.800    1
6.300 2.800 5.100 1.500    0
6.400 3.200 4.500 1.500    1
@EOF

chmod 640 NevProp/iris.net

exit 0
