PLearn 0.1
|
00001 // -*- C++ -*- 00002 00003 // CubicSpline.cc 00004 // 00005 // Copyright (C) 2007 Christian Dorion 00006 // 00007 // Redistribution and use in source and binary forms, with or without 00008 // modification, are permitted provided that the following conditions are met: 00009 // 00010 // 1. Redistributions of source code must retain the above copyright 00011 // notice, this list of conditions and the following disclaimer. 00012 // 00013 // 2. Redistributions in binary form must reproduce the above copyright 00014 // notice, this list of conditions and the following disclaimer in the 00015 // documentation and/or other materials provided with the distribution. 00016 // 00017 // 3. The name of the authors may not be used to endorse or promote 00018 // products derived from this software without specific prior written 00019 // permission. 00020 // 00021 // THIS SOFTWARE IS PROVIDED BY THE AUTHORS ``AS IS'' AND ANY EXPRESS OR 00022 // IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES 00023 // OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN 00024 // NO EVENT SHALL THE AUTHORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, 00025 // SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED 00026 // TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 00027 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 00028 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 00029 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 00030 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 00031 // 00032 // This file is part of the PLearn library. For more information on the PLearn 00033 // library, go to the PLearn Web site at www.plearn.org 00034 00035 // Authors: Christian Dorion 00036 00040 #include "CubicSpline.h" 00041 //#include "algo.h" 00042 00043 namespace PLearn { 00044 using namespace std; 00045 00046 PLEARN_IMPLEMENT_OBJECT( 00047 CubicSpline, 00048 "Unidimensional cubic spline learner.", 00049 "This learner fits a unidimensional cubic spline to a given set of\n" 00050 "points. That is, inputsize() must be one and the inputs are considered to be\n" 00051 "the x values, while targetsize() must also be one and the targets are\n" 00052 "considered to be the y values. The spline is fitted to the (x,y)-pairs so\n" 00053 "formed. X values don't need to be ordered; the ordering is ensured within\n" 00054 "the train method.\n"); 00055 00056 CubicSpline::CubicSpline() 00057 : m_low_slope(MISSING_VALUE), 00058 m_high_slope(MISSING_VALUE) 00059 { 00060 // Nothing to do here 00061 } 00062 00063 void CubicSpline::declareOptions(OptionList& ol) 00064 { 00065 declareOption(ol, "low_slope", &CubicSpline::m_low_slope, 00066 OptionBase::buildoption, 00067 "The slope to enforce at the leftmost node -- Default: NaN [None]"); 00068 00069 declareOption(ol, "high_slope", &CubicSpline::m_high_slope, 00070 OptionBase::buildoption, 00071 "The slope to enforce at the rightmost node -- Default: NaN [None]"); 00072 00073 // Learnt options 00074 declareOption(ol, "inputs", &CubicSpline::m_inputs, 00075 OptionBase::learntoption, 00076 "A buffer containing the last training set inputs"); 00077 00078 declareOption(ol, "targets", &CubicSpline::m_targets, 00079 OptionBase::learntoption, 00080 "A buffer containing the last training set targets"); 00081 00082 declareOption(ol, "coefficients", &CubicSpline::m_coefficients, 00083 OptionBase::learntoption, 00084 "The learnt coefficients"); 00085 00086 // Now call the parent class' declareOptions 00087 inherited::declareOptions(ol); 00088 } 00089 00090 void CubicSpline::build_() 00091 { 00092 // Nothing to do here 00093 } 00094 00095 // ### Nothing to add here, simply calls build_ 00096 void CubicSpline::build() 00097 { 00098 inherited::build(); 00099 build_(); 00100 } 00101 00102 00103 void CubicSpline::makeDeepCopyFromShallowCopy(CopiesMap& copies) 00104 { 00105 inherited::makeDeepCopyFromShallowCopy(copies); 00106 00107 deepCopyField(m_inputs, copies); 00108 deepCopyField(m_targets, copies); 00109 deepCopyField(m_coefficients, copies); 00110 } 00111 00112 00113 int CubicSpline::outputsize() const 00114 { 00115 return inputsize(); 00116 } 00117 00118 void CubicSpline::forget() 00119 { 00120 m_coefficients->resize(0); 00121 inherited::forget(); 00122 } 00123 00124 void CubicSpline::train() 00125 { 00126 // This learner fits unidimensional inputs/targets 00127 PLASSERT( inputsize() == 1 ); 00128 PLASSERT( targetsize() == 1 ); 00129 00130 // Train set is a member of PLearner; set through setTrainingSet() 00131 int n = train_set->length(); 00132 m_inputs = train_set.getColumn(0); 00133 m_targets = train_set.getColumn(1); 00134 PLASSERT( n >= 2 && train_set->width() == 2 ); 00135 00136 // Sort the inputs and targets along the inputs values 00137 TVec<int> indices = m_inputs.sortingPermutation(); 00138 m_inputs = m_inputs(indices); 00139 m_targets = m_targets(indices); 00140 00141 Vec u(n-1, 0.0); 00142 m_coefficients.resize(n); 00143 m_coefficients.fill(0.0); 00144 00145 // Low slope hack 00146 if ( !is_missing(m_low_slope) ) { 00147 u[0] = (3.0/(m_inputs[1]-m_inputs[0])) * 00148 ( (m_targets[1]-m_targets[0])/ 00149 (m_inputs[1] - m_inputs[0]) - m_low_slope ); 00150 m_coefficients[0] = -0.5; 00151 } 00152 00153 // Forward pass on coefficients 00154 for (int i=1; i < (n-1); i++) { 00155 real sig = (m_inputs[i]-m_inputs[i-1]) / (m_inputs[i+1]-m_inputs[i-1]); 00156 real p = sig*m_coefficients[i-1]+2.0; 00157 00158 m_coefficients[i] = (sig-1.0)/p; 00159 00160 u[i] = (m_targets[i+1]-m_targets[i]) / (m_inputs[i+1]-m_inputs[i]) 00161 - (m_targets[i]-m_targets[i-1]) / (m_inputs[i]-m_inputs[i-1]); 00162 00163 u[i] = (6.0*u[i]/(m_inputs[i+1]-m_inputs[i-1]) - sig*u[i-1]) / p; 00164 } 00165 00166 // High slope hack 00167 real un=0, qn=0; 00168 if ( !is_missing(m_high_slope) ) { 00169 qn = 0.5; 00170 un = (3.0/(m_inputs[n-1]-m_inputs[n-2])) * 00171 (m_high_slope - (m_targets[n-1]-m_targets[n-2]) / (m_inputs[n-1]-m_inputs[n-2])); 00172 } 00173 00174 // Compute the last coefficient 00175 m_coefficients[n-1] = (un-qn*u[n-2])/(qn*m_coefficients[n-1]+1.0); 00176 00177 // Backsubstitution step 00178 for (int k=n-2; k >= 0; k--) 00179 m_coefficients[k] = m_coefficients[k]*m_coefficients[k+1]+u[k]; 00180 } 00181 00182 00183 void CubicSpline::computeOutput(const Vec& input, Vec& output) const 00184 { 00185 PLASSERT( input.length() == 1 ); 00186 output.resize(1); 00187 00188 int n = m_inputs.length(); 00189 real x = input[0]; 00190 int pos = min( max(1, m_inputs.findSorted(x)), n-1 ); 00191 00192 real h = m_inputs[pos] - m_inputs[pos-1]; 00193 PLASSERT( h > 0.0 ); 00194 00195 real a = (m_inputs[pos] - x) / h; 00196 real b = (x - m_inputs[pos-1]) / h; 00197 output[0] = a*m_targets[pos-1] + b*m_targets[pos] 00198 + ((a*a*a - a)*m_coefficients[pos-1] 00199 + (b*b*b - b)*m_coefficients[pos]) * h*h/6.0; 00200 } 00201 00202 void CubicSpline::computeCostsFromOutputs(const Vec& input, const Vec& output, 00203 const Vec& target, Vec& costs) const 00204 { 00205 // No costs... 00206 } 00207 00208 TVec<string> CubicSpline::getTestCostNames() const 00209 { 00210 // None for now 00211 return TVec<string>(); 00212 } 00213 00214 TVec<string> CubicSpline::getTrainCostNames() const 00215 { 00216 // None for now 00217 return TVec<string>(); 00218 } 00219 00220 00221 } // end of namespace PLearn 00222 00223 00224 /* 00225 Local Variables: 00226 mode:c++ 00227 c-basic-offset:4 00228 c-file-style:"stroustrup" 00229 c-file-offsets:((innamespace . 0)(inline-open . 0)) 00230 indent-tabs-mode:nil 00231 fill-column:79 00232 End: 00233 */ 00234 // vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:encoding=utf-8:textwidth=79 :