Changeset 2415 for trunk/sources/LibSVM
- Timestamp:
- 10/07/09 11:58:21 (15 years ago)
- Location:
- trunk/sources/LibSVM
- Files:
-
- 3 added
- 1 deleted
- 17 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/sources/LibSVM/Cache.cs
r1819 r2415 49 49 private head_t lru_head; 50 50 51 internalCache(int count, long size)51 public Cache(int count, long size) 52 52 { 53 53 _count = count; … … 78 78 } 79 79 80 private static void swap<T>(ref T lhs, ref T rhs) 81 { 82 T tmp = lhs; 83 lhs = rhs; 84 rhs = tmp; 85 } 86 80 87 // request data [0,len) 81 88 // return some position p where [p,len) need to be filled 82 89 // (p >= len if nothing needs to be filled) 83 90 // java: simulate pointer using single-element array 84 internal virtual int get_data(int index, float[][] data, int len)91 public int GetData(int index, ref float[] data, int len) 85 92 { 86 93 head_t h = head[index]; … … 107 114 h.data = new_data; 108 115 _size -= more; 109 do 110 { 111 int _ = h.len; h.len = len; len = _; 112 } 113 while (false); 116 swap(ref h.len, ref len); 114 117 } 115 118 116 119 lru_insert(h); 117 data [0]= h.data;120 data = h.data; 118 121 return len; 119 122 } 120 123 121 internal virtual void swap_index(int i, int j)124 public void SwapIndex(int i, int j) 122 125 { 123 126 if (i == j) … … 128 131 if (head[j].len > 0) 129 132 lru_delete(head[j]); 130 do 131 { 132 float[] _ = head[i].data; head[i].data = head[j].data; head[j].data = _; 133 } 134 while (false); 135 do 136 { 137 int _ = head[i].len; head[i].len = head[j].len; head[j].len = _; 138 } 139 while (false); 133 swap(ref head[i].data, ref head[j].data); 134 swap(ref head[i].len, ref head[j].len); 140 135 if (head[i].len > 0) 141 136 lru_insert(head[i]); … … 144 139 145 140 if (i > j) 146 do 147 { 148 int _ = i; i = j; j = _; 149 } 150 while (false); 141 swap(ref i, ref j); 142 151 143 for (head_t h = lru_head.next; h != lru_head; h = h.next) 152 144 { … … 154 146 { 155 147 if (h.len > j) 156 do 157 { 158 float _ = h.data[i]; h.data[i] = h.data[j]; h.data[j] = _; 159 } 160 while (false); 148 swap(ref h.data[i], ref h.data[j]); 161 149 else 162 150 { -
trunk/sources/LibSVM/GaussianTransform.cs
r1819 r2415 21 21 using System.Collections.Generic; 22 22 using System.IO; 23 using System.Globalization; 24 using System.Threading; 23 25 24 26 namespace SVM 25 27 { 26 /// < remarks>28 /// <summary> 27 29 /// A transform which learns the mean and variance of a sample set and uses these to transform new data 28 30 /// so that it has zero mean and unit variance. 29 /// </ remarks>31 /// </summary> 30 32 public class GaussianTransform : IRangeTransform 31 33 { 32 private List<Node[]> _samples;33 private int _maxIndex;34 35 34 private double[] _means; 36 35 private double[] _stddevs; 37 36 38 37 /// <summary> 38 /// Determines the Gaussian transform for the provided problem. 39 /// </summary> 40 /// <param name="prob">The Problem to analyze</param> 41 /// <returns>The Gaussian transform for the problem</returns> 42 public static GaussianTransform Compute(Problem prob) 43 { 44 int[] counts = new int[prob.MaxIndex]; 45 double[] means = new double[prob.MaxIndex]; 46 foreach (Node[] sample in prob.X) 47 { 48 for (int i = 0; i < sample.Length; i++) 49 { 50 means[sample[i].Index-1] += sample[i].Value; 51 counts[sample[i].Index-1]++; 52 } 53 } 54 for (int i = 0; i < prob.MaxIndex; i++) 55 { 56 if (counts[i] == 0) 57 counts[i] = 2; 58 means[i] /= counts[i]; 59 } 60 61 double[] stddevs = new double[prob.MaxIndex]; 62 foreach (Node[] sample in prob.X) 63 { 64 for (int i = 0; i < sample.Length; i++) 65 { 66 double diff = sample[i].Value - means[sample[i].Index - 1]; 67 stddevs[sample[i].Index - 1] += diff * diff; 68 } 69 } 70 for (int i = 0; i < prob.MaxIndex; i++) 71 { 72 if (stddevs[i] == 0) 73 continue; 74 stddevs[i] /= (counts[i] - 1); 75 stddevs[i] = Math.Sqrt(stddevs[i]); 76 } 77 78 return new GaussianTransform(means, stddevs); 79 } 80 81 /// <summary> 39 82 /// Constructor. 40 83 /// </summary> 41 /// <param name="maxIndex">The maximum index of the vectors to be transformed</param> 42 public GaussianTransform(int maxIndex) 43 { 44 _samples = new List<Node[]>(); 45 } 46 private GaussianTransform(double[] means, double[] stddevs, int maxIndex) 84 /// <param name="means">Means in each dimension</param> 85 /// <param name="stddevs">Standard deviation in each dimension</param> 86 public GaussianTransform(double[] means, double[] stddevs) 47 87 { 48 88 _means = means; 49 89 _stddevs = stddevs; 50 _maxIndex = maxIndex;51 }52 53 /// <summary>54 /// Adds a sample to the data. No computation is performed. The maximum index of the55 /// sample must be less than MaxIndex.56 /// </summary>57 /// <param name="sample">The sample to add</param>58 public void Add(Node[] sample)59 {60 _samples.Add(sample);61 }62 63 /// <summary>64 /// Computes the statistics for the samples which have been obtained so far.65 /// </summary>66 public void ComputeStatistics()67 {68 int[] counts = new int[_maxIndex];69 _means = new double[_maxIndex];70 foreach(Node[] sample in _samples)71 {72 for (int i = 0; i < sample.Length; i++)73 {74 _means[sample[i].Index] += sample[i].Value;75 counts[sample[i].Index]++;76 }77 }78 for (int i = 0; i < _maxIndex; i++)79 {80 if (counts[i] == 0)81 counts[i] = 2;82 _means[i] /= counts[i];83 }84 85 _stddevs = new double[_maxIndex];86 foreach(Node[] sample in _samples)87 {88 for (int i = 0; i < sample.Length; i++)89 {90 double diff = sample[i].Value - _means[sample[i].Index];91 _stddevs[sample[i].Index] += diff * diff;92 }93 }94 for (int i = 0; i < _maxIndex; i++)95 {96 if (_stddevs[i] == 0)97 continue;98 _stddevs[i] /= (counts[i]-1);99 _stddevs[i] = Math.Sqrt(_stddevs[i]);100 }101 90 } 102 91 … … 109 98 public static void Write(Stream stream, GaussianTransform transform) 110 99 { 100 TemporaryCulture.Start(); 101 111 102 StreamWriter output = new StreamWriter(stream); 112 output.WriteLine(transform._m axIndex);113 for (int i = 0; i < transform._m axIndex; i++)103 output.WriteLine(transform._means.Length); 104 for (int i = 0; i < transform._means.Length; i++) 114 105 output.WriteLine("{0} {1}", transform._means[i], transform._stddevs[i]); 115 106 output.Flush(); 107 108 TemporaryCulture.Stop(); 116 109 } 117 110 … … 123 116 public static GaussianTransform Read(Stream stream) 124 117 { 118 TemporaryCulture.Start(); 119 125 120 StreamReader input = new StreamReader(stream); 126 int length = int.Parse(input.ReadLine() );121 int length = int.Parse(input.ReadLine(), CultureInfo.InvariantCulture); 127 122 double[] means = new double[length]; 128 123 double[] stddevs = new double[length]; … … 130 125 { 131 126 string[] parts = input.ReadLine().Split(); 132 means[i] = double.Parse(parts[0]); 133 stddevs[i] = double.Parse(parts[1]); 134 } 135 return new GaussianTransform(means, stddevs, length); 127 means[i] = double.Parse(parts[0], CultureInfo.InvariantCulture); 128 stddevs[i] = double.Parse(parts[1], CultureInfo.InvariantCulture); 129 } 130 131 TemporaryCulture.Stop(); 132 133 return new GaussianTransform(means, stddevs); 136 134 } 137 135 … … 177 175 /// <summary> 178 176 /// Transform the input value using the transform stored for the provided index. 179 /// <see cref="ComputeStatistics"/> must be called first, or the transform must180 /// have been read from the disk.181 177 /// </summary> 182 178 /// <param name="input">Input value</param> … … 185 181 public double Transform(double input, int index) 186 182 { 183 index--; 187 184 if (_stddevs[index] == 0) 188 185 return 0; … … 192 189 } 193 190 /// <summary> 194 /// Transforms the input array. <see cref="ComputeStatistics"/> must be called 195 /// first, or the transform must have been read from the disk. 191 /// Transforms the input array. 196 192 /// </summary> 197 193 /// <param name="input">The array to transform</param> -
trunk/sources/LibSVM/IRangeTransform.cs
r1819 r2415 20 20 namespace SVM 21 21 { 22 /// < remarks>22 /// <summary> 23 23 /// Interface implemented by range transforms. 24 /// </ remarks>24 /// </summary> 25 25 public interface IRangeTransform 26 26 { -
trunk/sources/LibSVM/LibSVM.csproj
r1846 r2415 80 80 </ItemGroup> 81 81 <ItemGroup> 82 <Compile Include="Kernel.cs" /> 82 83 <Compile Include="Properties\AssemblyInfo.cs" /> 83 84 <Compile Include="Cache.cs" /> … … 95 96 <Compile Include="Scaling.cs" /> 96 97 <Compile Include="Solver.cs" /> 97 <Compile Include="SupportClass.cs" /> 98 <Compile Include="SVMExtensions.cs" /> 99 <Compile Include="TemporaryCulture.cs" /> 98 100 <Compile Include="Training.cs" /> 99 101 </ItemGroup> -
trunk/sources/LibSVM/Model.cs
r2411 r2415 21 21 using System; 22 22 using System.IO; 23 24 namespace SVM { 25 /// <remarks> 26 /// Encapsulates an SVM Model. 27 /// </remarks> 28 [Serializable] 29 public class Model { 30 private Parameter _parameter; 31 private int _numberOfClasses; 32 private int _supportVectorCount; 33 private Node[][] _supportVectors; 34 private double[][] _supportVectorCoefficients; 35 private double[] _rho; 36 private double[] _pairwiseProbabilityA; 37 private double[] _pairwiseProbabilityB; 38 39 private int[] _classLabels; 40 private int[] _numberOfSVPerClass; 41 42 internal Model() { 43 } 44 23 using System.Threading; 24 using System.Globalization; 25 26 namespace SVM 27 { 45 28 /// <summary> 46 /// Parameter object.29 /// Encapsulates an SVM Model. 47 30 /// </summary> 48 public Parameter Parameter { 49 get { 50 return _parameter; 51 } 52 set { 53 _parameter = value; 54 } 55 } 56 57 /// <summary> 58 /// Number of classes in the model. 59 /// </summary> 60 public int NumberOfClasses { 61 get { 62 return _numberOfClasses; 63 } 64 set { 65 _numberOfClasses = value; 66 } 67 } 68 69 /// <summary> 70 /// Total number of support vectors. 71 /// </summary> 72 public int SupportVectorCount { 73 get { 74 return _supportVectorCount; 75 } 76 set { 77 _supportVectorCount = value; 78 } 79 } 80 81 /// <summary> 82 /// The support vectors. 83 /// </summary> 84 public Node[][] SupportVectors { 85 get { 86 return _supportVectors; 87 } 88 set { 89 _supportVectors = value; 90 } 91 } 92 93 /// <summary> 94 /// The coefficients for the support vectors. 95 /// </summary> 96 public double[][] SupportVectorCoefficients { 97 get { 98 return _supportVectorCoefficients; 99 } 100 set { 101 _supportVectorCoefficients = value; 102 } 103 } 104 105 /// <summary> 106 /// Rho values. 107 /// </summary> 108 public double[] Rho { 109 get { 110 return _rho; 111 } 112 set { 113 _rho = value; 114 } 115 } 116 117 /// <summary> 118 /// First pairwise probability. 119 /// </summary> 120 public double[] PairwiseProbabilityA { 121 get { 122 return _pairwiseProbabilityA; 123 } 124 set { 125 _pairwiseProbabilityA = value; 126 } 127 } 128 129 /// <summary> 130 /// Second pairwise probability. 131 /// </summary> 132 public double[] PairwiseProbabilityB { 133 get { 134 return _pairwiseProbabilityB; 135 } 136 set { 137 _pairwiseProbabilityB = value; 138 } 139 } 140 141 // for classification only 142 143 /// <summary> 144 /// Class labels. 145 /// </summary> 146 public int[] ClassLabels { 147 get { 148 return _classLabels; 149 } 150 set { 151 _classLabels = value; 152 } 153 } 154 155 /// <summary> 156 /// Number of support vectors per class. 157 /// </summary> 158 public int[] NumberOfSVPerClass { 159 get { 160 return _numberOfSVPerClass; 161 } 162 set { 163 _numberOfSVPerClass = value; 164 } 165 } 166 167 /// <summary> 168 /// Reads a Model from the provided file. 169 /// </summary> 170 /// <param name="filename">The name of the file containing the Model</param> 171 /// <returns>the Model</returns> 172 public static Model Read(string filename) { 173 FileStream input = File.OpenRead(filename); 174 try { 175 return Read(input); 176 } 177 finally { 178 input.Close(); 179 } 180 } 181 182 /// <summary> 183 /// Reads a Model from the provided stream. 184 /// </summary> 185 /// <param name="stream">The stream from which to read the Model.</param> 186 /// <returns>the Model</returns> 187 public static Model Read(Stream stream) { 188 StreamReader input = new StreamReader(stream); 189 190 // read parameters 191 192 Model model = new Model(); 193 Parameter param = new Parameter(); 194 model.Parameter = param; 195 model.Rho = null; 196 model.PairwiseProbabilityA = null; 197 model.PairwiseProbabilityB = null; 198 model.ClassLabels = null; 199 model.NumberOfSVPerClass = null; 200 201 bool headerFinished = false; 202 while (!headerFinished) { 203 string line = input.ReadLine(); 204 string cmd, arg; 205 int splitIndex = line.IndexOf(' '); 206 if (splitIndex >= 0) { 207 cmd = line.Substring(0, splitIndex); 208 arg = line.Substring(splitIndex + 1); 209 } else { 210 cmd = line; 211 arg = ""; 212 } 213 arg = arg.ToLower(); 214 215 int i, n; 216 switch (cmd) { 217 case "svm_type": 218 param.SvmType = (SvmType)Enum.Parse(typeof(SvmType), arg.ToUpper()); 219 break; 220 221 case "kernel_type": 222 param.KernelType = (KernelType)Enum.Parse(typeof(KernelType), arg.ToUpper()); 223 break; 224 225 case "degree": 226 param.Degree = int.Parse(arg); 227 break; 228 229 case "gamma": 230 param.Gamma = double.Parse(arg); 231 break; 232 233 case "coef0": 234 param.Coefficient0 = double.Parse(arg); 235 break; 236 237 case "nr_class": 238 model.NumberOfClasses = int.Parse(arg); 239 break; 240 241 case "total_sv": 242 model.SupportVectorCount = int.Parse(arg); 243 break; 244 245 case "rho": 246 n = model.NumberOfClasses * (model.NumberOfClasses - 1) / 2; 247 model.Rho = new double[n]; 248 string[] rhoParts = arg.Split(); 249 for (i = 0; i < n; i++) 250 model.Rho[i] = double.Parse(rhoParts[i]); 251 break; 252 253 case "label": 254 n = model.NumberOfClasses; 255 model.ClassLabels = new int[n]; 256 string[] labelParts = arg.Split(); 257 for (i = 0; i < n; i++) 258 model.ClassLabels[i] = int.Parse(labelParts[i]); 259 break; 260 261 case "probA": 262 n = model.NumberOfClasses * (model.NumberOfClasses - 1) / 2; 263 model.PairwiseProbabilityA = new double[n]; 264 string[] probAParts = arg.Split(); 265 for (i = 0; i < n; i++) 266 model.PairwiseProbabilityA[i] = double.Parse(probAParts[i]); 267 break; 268 269 case "probB": 270 n = model.NumberOfClasses * (model.NumberOfClasses - 1) / 2; 271 model.PairwiseProbabilityB = new double[n]; 272 string[] probBParts = arg.Split(); 273 for (i = 0; i < n; i++) 274 model.PairwiseProbabilityB[i] = double.Parse(probBParts[i]); 275 break; 276 277 case "nr_sv": 278 n = model.NumberOfClasses; 279 model.NumberOfSVPerClass = new int[n]; 280 string[] nrsvParts = arg.Split(); 281 for (i = 0; i < n; i++) 282 model.NumberOfSVPerClass[i] = int.Parse(nrsvParts[i]); 283 break; 284 285 case "SV": 286 headerFinished = true; 287 break; 288 289 default: 290 throw new Exception("Unknown text in model file"); 291 } 292 } 293 294 // read sv_coef and SV 295 296 int m = model.NumberOfClasses - 1; 297 int l = model.SupportVectorCount; 298 model.SupportVectorCoefficients = new double[m][]; 299 for (int i = 0; i < m; i++) { 300 model.SupportVectorCoefficients[i] = new double[l]; 301 } 302 model.SupportVectors = new Node[l][]; 303 304 for (int i = 0; i < l; i++) { 305 string[] parts = input.ReadLine().Trim().Split(); 306 307 for (int k = 0; k < m; k++) 308 model.SupportVectorCoefficients[k][i] = double.Parse(parts[k]); 309 int n = parts.Length - m; 310 model.SupportVectors[i] = new Node[n]; 311 for (int j = 0; j < n; j++) { 312 string[] nodeParts = parts[m + j].Split(':'); 313 model.SupportVectors[i][j] = new Node(); 314 model.SupportVectors[i][j].Index = int.Parse(nodeParts[0]); 315 model.SupportVectors[i][j].Value = double.Parse(nodeParts[1]); 316 } 317 } 318 319 return model; 320 } 321 322 /// <summary> 323 /// Writes a model to the provided filename. This will overwrite any previous data in the file. 324 /// </summary> 325 /// <param name="filename">The desired file</param> 326 /// <param name="model">The Model to write</param> 327 public static void Write(string filename, Model model) { 328 FileStream stream = File.Open(filename, FileMode.Create); 329 try { 330 Write(stream, model); 331 } 332 finally { 333 stream.Close(); 334 } 335 } 336 337 /// <summary> 338 /// Writes a model to the provided stream. 339 /// </summary> 340 /// <param name="stream">The output stream</param> 341 /// <param name="model">The model to write</param> 342 public static void Write(Stream stream, Model model) { 343 StreamWriter output = new StreamWriter(stream); 344 345 Parameter param = model.Parameter; 346 347 output.Write("svm_type " + param.SvmType + Environment.NewLine); 348 output.Write("kernel_type " + param.KernelType + Environment.NewLine); 349 350 if (param.KernelType == KernelType.POLY) 351 output.Write("degree " + param.Degree + Environment.NewLine); 352 353 if (param.KernelType == KernelType.POLY || param.KernelType == KernelType.RBF || param.KernelType == KernelType.SIGMOID) 354 output.Write("gamma " + param.Gamma + Environment.NewLine); 355 356 if (param.KernelType == KernelType.POLY || param.KernelType == KernelType.SIGMOID) 357 output.Write("coef0 " + param.Coefficient0 + Environment.NewLine); 358 359 int nr_class = model.NumberOfClasses; 360 int l = model.SupportVectorCount; 361 output.Write("nr_class " + nr_class + Environment.NewLine); 362 output.Write("total_sv " + l + Environment.NewLine); 363 364 { 365 output.Write("rho"); 366 for (int i = 0; i < nr_class * (nr_class - 1) / 2; i++) 367 output.Write(" " + model.Rho[i]); 368 output.Write(Environment.NewLine); 369 } 370 371 if (model.ClassLabels != null) { 372 output.Write("label"); 373 for (int i = 0; i < nr_class; i++) 374 output.Write(" " + model.ClassLabels[i]); 375 output.Write(Environment.NewLine); 376 } 377 378 if (model.PairwiseProbabilityA != null) 379 // regression has probA only 380 { 381 output.Write("probA"); 382 for (int i = 0; i < nr_class * (nr_class - 1) / 2; i++) 383 output.Write(" " + model.PairwiseProbabilityA[i]); 384 output.Write(Environment.NewLine); 385 } 386 if (model.PairwiseProbabilityB != null) { 387 output.Write("probB"); 388 for (int i = 0; i < nr_class * (nr_class - 1) / 2; i++) 389 output.Write(" " + model.PairwiseProbabilityB[i]); 390 output.Write(Environment.NewLine); 391 } 392 393 if (model.NumberOfSVPerClass != null) { 394 output.Write("nr_sv"); 395 for (int i = 0; i < nr_class; i++) 396 output.Write(" " + model.NumberOfSVPerClass[i]); 397 output.Write(Environment.NewLine); 398 } 399 400 output.Write("SV\n"); 401 double[][] sv_coef = model.SupportVectorCoefficients; 402 Node[][] SV = model.SupportVectors; 403 404 for (int i = 0; i < l; i++) { 405 for (int j = 0; j < nr_class - 1; j++) 406 output.Write(sv_coef[j][i] + " "); 407 408 Node[] p = SV[i]; 409 if (p.Length == 0) { 410 output.WriteLine(); 411 continue; 412 } 413 if (param.KernelType == KernelType.PRECOMPUTED) 414 output.Write("0:{0}", (int)p[0].Value); 415 else { 416 output.Write("{0}:{1}", p[0].Index, p[0].Value); 417 for (int j = 1; j < p.Length; j++) 418 output.Write(" {0}:{1}", p[j].Index, p[j].Value); 419 } 420 output.WriteLine(); 421 } 422 423 output.Flush(); 424 } 425 } 31 [Serializable] 32 public class Model 33 { 34 private Parameter _parameter; 35 private int _numberOfClasses; 36 private int _supportVectorCount; 37 private Node[][] _supportVectors; 38 private double[][] _supportVectorCoefficients; 39 private double[] _rho; 40 private double[] _pairwiseProbabilityA; 41 private double[] _pairwiseProbabilityB; 42 43 private int[] _classLabels; 44 private int[] _numberOfSVPerClass; 45 46 internal Model() 47 { 48 } 49 50 /// <summary> 51 /// Parameter object. 52 /// </summary> 53 public Parameter Parameter 54 { 55 get 56 { 57 return _parameter; 58 } 59 set 60 { 61 _parameter = value; 62 } 63 } 64 65 /// <summary> 66 /// Number of classes in the model. 67 /// </summary> 68 public int NumberOfClasses 69 { 70 get 71 { 72 return _numberOfClasses; 73 } 74 set 75 { 76 _numberOfClasses = value; 77 } 78 } 79 80 /// <summary> 81 /// Total number of support vectors. 82 /// </summary> 83 public int SupportVectorCount 84 { 85 get 86 { 87 return _supportVectorCount; 88 } 89 set 90 { 91 _supportVectorCount = value; 92 } 93 } 94 95 /// <summary> 96 /// The support vectors. 97 /// </summary> 98 public Node[][] SupportVectors 99 { 100 get 101 { 102 return _supportVectors; 103 } 104 set 105 { 106 _supportVectors = value; 107 } 108 } 109 110 /// <summary> 111 /// The coefficients for the support vectors. 112 /// </summary> 113 public double[][] SupportVectorCoefficients 114 { 115 get 116 { 117 return _supportVectorCoefficients; 118 } 119 set 120 { 121 _supportVectorCoefficients = value; 122 } 123 } 124 125 /// <summary> 126 /// Rho values. 127 /// </summary> 128 public double[] Rho 129 { 130 get 131 { 132 return _rho; 133 } 134 set 135 { 136 _rho = value; 137 } 138 } 139 140 /// <summary> 141 /// First pairwise probability. 142 /// </summary> 143 public double[] PairwiseProbabilityA 144 { 145 get 146 { 147 return _pairwiseProbabilityA; 148 } 149 set 150 { 151 _pairwiseProbabilityA = value; 152 } 153 } 154 155 /// <summary> 156 /// Second pairwise probability. 157 /// </summary> 158 public double[] PairwiseProbabilityB 159 { 160 get 161 { 162 return _pairwiseProbabilityB; 163 } 164 set 165 { 166 _pairwiseProbabilityB = value; 167 } 168 } 169 170 // for classification only 171 172 /// <summary> 173 /// Class labels. 174 /// </summary> 175 public int[] ClassLabels 176 { 177 get 178 { 179 return _classLabels; 180 } 181 set 182 { 183 _classLabels = value; 184 } 185 } 186 187 /// <summary> 188 /// Number of support vectors per class. 189 /// </summary> 190 public int[] NumberOfSVPerClass 191 { 192 get 193 { 194 return _numberOfSVPerClass; 195 } 196 set 197 { 198 _numberOfSVPerClass = value; 199 } 200 } 201 202 /// <summary> 203 /// Reads a Model from the provided file. 204 /// </summary> 205 /// <param name="filename">The name of the file containing the Model</param> 206 /// <returns>the Model</returns> 207 public static Model Read(string filename) 208 { 209 FileStream input = File.OpenRead(filename); 210 try 211 { 212 return Read(input); 213 } 214 finally 215 { 216 input.Close(); 217 } 218 } 219 220 /// <summary> 221 /// Reads a Model from the provided stream. 222 /// </summary> 223 /// <param name="stream">The stream from which to read the Model.</param> 224 /// <returns>the Model</returns> 225 public static Model Read(Stream stream) 226 { 227 TemporaryCulture.Start(); 228 229 StreamReader input = new StreamReader(stream); 230 231 // read parameters 232 233 Model model = new Model(); 234 Parameter param = new Parameter(); 235 model.Parameter = param; 236 model.Rho = null; 237 model.PairwiseProbabilityA = null; 238 model.PairwiseProbabilityB = null; 239 model.ClassLabels = null; 240 model.NumberOfSVPerClass = null; 241 242 bool headerFinished = false; 243 while (!headerFinished) 244 { 245 string line = input.ReadLine(); 246 string cmd, arg; 247 int splitIndex = line.IndexOf(' '); 248 if (splitIndex >= 0) 249 { 250 cmd = line.Substring(0, splitIndex); 251 arg = line.Substring(splitIndex + 1); 252 } 253 else 254 { 255 cmd = line; 256 arg = ""; 257 } 258 arg = arg.ToLower(); 259 260 int i,n; 261 switch(cmd){ 262 case "svm_type": 263 param.SvmType = (SvmType)Enum.Parse(typeof(SvmType), arg.ToUpper()); 264 break; 265 266 case "kernel_type": 267 param.KernelType = (KernelType)Enum.Parse(typeof(KernelType), arg.ToUpper()); 268 break; 269 270 case "degree": 271 param.Degree = int.Parse(arg); 272 break; 273 274 case "gamma": 275 param.Gamma = double.Parse(arg); 276 break; 277 278 case "coef0": 279 param.Coefficient0 = double.Parse(arg); 280 break; 281 282 case "nr_class": 283 model.NumberOfClasses = int.Parse(arg); 284 break; 285 286 case "total_sv": 287 model.SupportVectorCount = int.Parse(arg); 288 break; 289 290 case "rho": 291 n = model.NumberOfClasses * (model.NumberOfClasses - 1) / 2; 292 model.Rho = new double[n]; 293 string[] rhoParts = arg.Split(); 294 for(i=0; i<n; i++) 295 model.Rho[i] = double.Parse(rhoParts[i]); 296 break; 297 298 case "label": 299 n = model.NumberOfClasses; 300 model.ClassLabels = new int[n]; 301 string[] labelParts = arg.Split(); 302 for (i = 0; i < n; i++) 303 model.ClassLabels[i] = int.Parse(labelParts[i]); 304 break; 305 306 case "probA": 307 n = model.NumberOfClasses * (model.NumberOfClasses - 1) / 2; 308 model.PairwiseProbabilityA = new double[n]; 309 string[] probAParts = arg.Split(); 310 for (i = 0; i < n; i++) 311 model.PairwiseProbabilityA[i] = double.Parse(probAParts[i]); 312 break; 313 314 case "probB": 315 n = model.NumberOfClasses * (model.NumberOfClasses - 1) / 2; 316 model.PairwiseProbabilityB = new double[n]; 317 string[] probBParts = arg.Split(); 318 for (i = 0; i < n; i++) 319 model.PairwiseProbabilityB[i] = double.Parse(probBParts[i]); 320 break; 321 322 case "nr_sv": 323 n = model.NumberOfClasses; 324 model.NumberOfSVPerClass = new int[n]; 325 string[] nrsvParts = arg.Split(); 326 for (i = 0; i < n; i++) 327 model.NumberOfSVPerClass[i] = int.Parse(nrsvParts[i]); 328 break; 329 330 case "SV": 331 headerFinished = true; 332 break; 333 334 default: 335 throw new Exception("Unknown text in model file"); 336 } 337 } 338 339 // read sv_coef and SV 340 341 int m = model.NumberOfClasses - 1; 342 int l = model.SupportVectorCount; 343 model.SupportVectorCoefficients = new double[m][]; 344 for (int i = 0; i < m; i++) 345 { 346 model.SupportVectorCoefficients[i] = new double[l]; 347 } 348 model.SupportVectors = new Node[l][]; 349 350 for (int i = 0; i < l; i++) 351 { 352 string[] parts = input.ReadLine().Trim().Split(); 353 354 for (int k = 0; k < m; k++) 355 model.SupportVectorCoefficients[k][i] = double.Parse(parts[k]); 356 int n = parts.Length-m; 357 model.SupportVectors[i] = new Node[n]; 358 for (int j = 0; j < n; j++) 359 { 360 string[] nodeParts = parts[m + j].Split(':'); 361 model.SupportVectors[i][j] = new Node(); 362 model.SupportVectors[i][j].Index = int.Parse(nodeParts[0]); 363 model.SupportVectors[i][j].Value = double.Parse(nodeParts[1]); 364 } 365 } 366 367 TemporaryCulture.Stop(); 368 369 return model; 370 } 371 372 /// <summary> 373 /// Writes a model to the provided filename. This will overwrite any previous data in the file. 374 /// </summary> 375 /// <param name="filename">The desired file</param> 376 /// <param name="model">The Model to write</param> 377 public static void Write(string filename, Model model) 378 { 379 FileStream stream = File.Open(filename, FileMode.Create); 380 try 381 { 382 Write(stream, model); 383 } 384 finally 385 { 386 stream.Close(); 387 } 388 } 389 390 /// <summary> 391 /// Writes a model to the provided stream. 392 /// </summary> 393 /// <param name="stream">The output stream</param> 394 /// <param name="model">The model to write</param> 395 public static void Write(Stream stream, Model model) 396 { 397 TemporaryCulture.Start(); 398 399 StreamWriter output = new StreamWriter(stream); 400 401 Parameter param = model.Parameter; 402 403 output.Write("svm_type " + param.SvmType + "\n"); 404 output.Write("kernel_type " + param.KernelType + "\n"); 405 406 if (param.KernelType == KernelType.POLY) 407 output.Write("degree " + param.Degree + "\n"); 408 409 if (param.KernelType == KernelType.POLY || param.KernelType == KernelType.RBF || param.KernelType == KernelType.SIGMOID) 410 output.Write("gamma " + param.Gamma + "\n"); 411 412 if (param.KernelType == KernelType.POLY || param.KernelType == KernelType.SIGMOID) 413 output.Write("coef0 " + param.Coefficient0 + "\n"); 414 415 int nr_class = model.NumberOfClasses; 416 int l = model.SupportVectorCount; 417 output.Write("nr_class " + nr_class + "\n"); 418 output.Write("total_sv " + l + "\n"); 419 420 { 421 output.Write("rho"); 422 for (int i = 0; i < nr_class * (nr_class - 1) / 2; i++) 423 output.Write(" " + model.Rho[i]); 424 output.Write("\n"); 425 } 426 427 if (model.ClassLabels != null) 428 { 429 output.Write("label"); 430 for (int i = 0; i < nr_class; i++) 431 output.Write(" " + model.ClassLabels[i]); 432 output.Write("\n"); 433 } 434 435 if (model.PairwiseProbabilityA != null) 436 // regression has probA only 437 { 438 output.Write("probA"); 439 for (int i = 0; i < nr_class * (nr_class - 1) / 2; i++) 440 output.Write(" " + model.PairwiseProbabilityA[i]); 441 output.Write("\n"); 442 } 443 if (model.PairwiseProbabilityB != null) 444 { 445 output.Write("probB"); 446 for (int i = 0; i < nr_class * (nr_class - 1) / 2; i++) 447 output.Write(" " + model.PairwiseProbabilityB[i]); 448 output.Write("\n"); 449 } 450 451 if (model.NumberOfSVPerClass != null) 452 { 453 output.Write("nr_sv"); 454 for (int i = 0; i < nr_class; i++) 455 output.Write(" " + model.NumberOfSVPerClass[i]); 456 output.Write("\n"); 457 } 458 459 output.Write("SV\n"); 460 double[][] sv_coef = model.SupportVectorCoefficients; 461 Node[][] SV = model.SupportVectors; 462 463 for (int i = 0; i < l; i++) 464 { 465 for (int j = 0; j < nr_class - 1; j++) 466 output.Write(sv_coef[j][i] + " "); 467 468 Node[] p = SV[i]; 469 if (p.Length == 0) 470 { 471 output.WriteLine(); 472 continue; 473 } 474 if (param.KernelType == KernelType.PRECOMPUTED) 475 output.Write("0:{0}", (int)p[0].Value); 476 else 477 { 478 output.Write("{0}:{1}", p[0].Index, p[0].Value); 479 for (int j = 1; j < p.Length; j++) 480 output.Write(" {0}:{1}", p[j].Index, p[j].Value); 481 } 482 output.WriteLine(); 483 } 484 485 output.Flush(); 486 487 TemporaryCulture.Stop(); 488 } 489 } 426 490 } -
trunk/sources/LibSVM/Node.cs
r1819 r2415 22 22 namespace SVM 23 23 { 24 /// < remarks>24 /// <summary> 25 25 /// Encapsulates a node in a Problem vector, with an index and a value (for more efficient representation 26 26 /// of sparse data. 27 /// </ remarks>27 /// </summary> 28 28 [Serializable] 29 29 public class Node : IComparable<Node> 30 30 { 31 privateint _index;32 privatedouble _value;31 internal int _index; 32 internal double _value; 33 33 34 34 /// <summary> -
trunk/sources/LibSVM/Parameter.cs
r1819 r2415 19 19 20 20 using System; 21 using System.Linq; 22 using System.Collections.Generic; 21 23 22 24 namespace SVM 23 25 { 24 /// < remarks>26 /// <summary> 25 27 /// Contains all of the types of SVM this library can model. 26 /// </ remarks>28 /// </summary> 27 29 public enum SvmType { 28 30 /// <summary> … … 47 49 NU_SVR 48 50 }; 49 /// < remarks>51 /// <summary> 50 52 /// Contains the various kernel types this library can use. 51 /// </ remarks>53 /// </summary> 52 54 public enum KernelType { 53 55 /// <summary> … … 73 75 }; 74 76 75 /// < remarks>77 /// <summary> 76 78 /// This class contains the various parameters which can affect the way in which an SVM 77 79 /// is learned. Unless you know what you are doing, chances are you are best off using 78 80 /// the default values. 79 /// </ remarks>81 /// </summary> 80 82 [Serializable] 81 83 public class Parameter : ICloneable … … 91 93 private double _eps; 92 94 93 private int _weightCount; 94 private int[] _weightLabels; 95 private double[] _weights; 95 private Dictionary<int, double> _weights; 96 96 private double _nu; 97 97 private double _p; … … 116 116 _shrinking = true; 117 117 _probability = false; 118 _weightCount = 0; 119 _weightLabels = new int[0]; 120 _weights = new double[0]; 118 _weights = new Dictionary<int, double>(); 121 119 } 122 120 … … 234 232 } 235 233 } 236 /// <summary> 237 /// Number of weights. 238 /// </summary> 239 public int WeightCount 240 { 241 get 242 { 243 return _weightCount; 244 } 245 set 246 { 247 _weightCount = value; 248 } 249 } 250 /// <summary> 251 /// Array of indicies corresponding to the Weights array (for C-SVC) 252 /// </summary> 253 public int[] WeightLabels 254 { 255 get 256 { 257 return _weightLabels; 258 } 259 set 260 { 261 _weightLabels = value; 262 } 263 } 264 /// <summary> 265 /// The parameter C of class i to weight*C in C-SVC (default 1) 266 /// </summary> 267 public double[] Weights 268 { 269 get 270 { 234 235 /// <summary> 236 /// Contains custom weights for class labels. Default weight value is 1. 237 /// </summary> 238 public Dictionary<int,double> Weights 239 { 240 get{ 271 241 return _weights; 272 242 } 273 set 274 { 275 _weights = value; 276 } 277 } 243 } 244 278 245 /// <summary> 279 246 /// The parameter nu of nu-SVC, one-class SVM, and nu-SVR (default 0.5) -
trunk/sources/LibSVM/ParameterSelection.cs
r1819 r2415 24 24 namespace SVM 25 25 { 26 /// < remarks>26 /// <summary> 27 27 /// This class contains routines which perform parameter selection for a model which uses C-SVC and 28 28 /// an RBF kernel. 29 /// </ remarks>29 /// </summary> 30 30 public static class ParameterSelection 31 31 { … … 142 142 Gamma = 0; 143 143 double crossValidation = double.MinValue; 144 StreamWriter output = new StreamWriter("graph.txt"); 144 StreamWriter output = null; 145 if(outputFile != null) 146 output = new StreamWriter(outputFile); 145 147 for(int i=0; i<CValues.Count; i++) 146 148 for (int j = 0; j < GammaValues.Count; j++) … … 150 152 double test = Training.PerformCrossValidation(problem, parameters, nrfold); 151 153 Console.Write("{0} {1} {2}", parameters.C, parameters.Gamma, test); 152 output.WriteLine("{0} {1} {2}", parameters.C, parameters.Gamma, test); 154 if(output != null) 155 output.WriteLine("{0} {1} {2}", parameters.C, parameters.Gamma, test); 153 156 if (test > crossValidation) 154 157 { … … 160 163 else Console.WriteLine(); 161 164 } 162 output.Close(); 165 if(output != null) 166 output.Close(); 163 167 } 164 168 /// <summary> … … 207 211 Gamma = 0; 208 212 double maxScore = double.MinValue; 209 StreamWriter output = new StreamWriter(outputFile); 213 StreamWriter output = null; 214 if(outputFile != null) 215 output = new StreamWriter(outputFile); 210 216 for (int i = 0; i < CValues.Count; i++) 211 217 for (int j = 0; j < GammaValues.Count; j++) … … 216 222 double test = Prediction.Predict(validation, "tmp.txt", model, false); 217 223 Console.Write("{0} {1} {2}", parameters.C, parameters.Gamma, test); 218 output.WriteLine("{0} {1} {2}", parameters.C, parameters.Gamma, test); 224 if(output != null) 225 output.WriteLine("{0} {1} {2}", parameters.C, parameters.Gamma, test); 219 226 if (test > maxScore) 220 227 { … … 226 233 else Console.WriteLine(); 227 234 } 228 output.Close(); 235 if(output != null) 236 output.Close(); 229 237 } 230 238 } -
trunk/sources/LibSVM/PerformanceEvaluator.cs
r1819 r2415 21 21 using System.Collections.Generic; 22 22 using System.IO; 23 using System.Globalization; 23 24 24 25 namespace SVM 25 26 { 26 /// < remarks>27 /// <summary> 27 28 /// Class encoding a member of a ranked set of labels. 28 /// </ remarks>29 /// </summary> 29 30 public class RankPair : IComparable<RankPair> 30 31 { … … 138 139 } 139 140 140 /// < remarks>141 /// <summary> 141 142 /// Class which evaluates an SVM model using several standard techniques. 142 /// </ remarks>143 /// </summary> 143 144 public class PerformanceEvaluator 144 145 { … … 185 186 /// <param name="model">Model to evaluate</param> 186 187 /// <param name="problem">Problem to evaluate</param> 187 /// <param name=" label">Label to be evaluate for</param>188 public PerformanceEvaluator(Model model, Problem problem, double label) : this(model, problem, label, "tmp.results") { }188 /// <param name="category">Label to be evaluate for</param> 189 public PerformanceEvaluator(Model model, Problem problem, double category) : this(model, problem, category, "tmp.results") { } 189 190 /// <summary> 190 191 /// Constructor. … … 220 221 int confidenceIndex = -1; 221 222 for (int i = 1; i < parts.Length; i++) 222 if (double.Parse(parts[i] ) == category)223 if (double.Parse(parts[i], CultureInfo.InvariantCulture) == category) 223 224 { 224 225 confidenceIndex = i; … … 229 230 { 230 231 parts = input.ReadLine().Split(new char[] { ' ' }, StringSplitOptions.RemoveEmptyEntries); 231 double confidence = double.Parse(parts[confidenceIndex] );232 _data.Add(new RankPair(confidence, labels[i] ));232 double confidence = double.Parse(parts[confidenceIndex], CultureInfo.InvariantCulture); 233 _data.Add(new RankPair(confidence, labels[i] == category ? 1 : 0)); 233 234 } 234 235 input.Close(); -
trunk/sources/LibSVM/PrecomputedKernel.cs
r1819 r2415 23 23 namespace SVM 24 24 { 25 /// < remarks>25 /// <summary> 26 26 /// Class encapsulating a precomputed kernel, where each position indicates the similarity score for two items in the training data. 27 /// </ remarks>27 /// </summary> 28 28 [Serializable] 29 29 public class PrecomputedKernel … … 42 42 _rows = _similarities.GetLength(0); 43 43 _columns = _similarities.GetLength(1); 44 } 45 46 /// <summary> 47 /// Constructor. 48 /// </summary> 49 /// <param name="nodes">Nodes for self-similarity analysis</param> 50 /// <param name="param">Parameters to use when computing similarities</param> 51 public PrecomputedKernel(List<Node[]> nodes, Parameter param) 52 { 53 _rows = nodes.Count; 54 _columns = _rows; 55 _similarities = new float[_rows, _columns]; 56 for (int r = 0; r < _rows; r++) 57 { 58 for (int c = 0; c < r; c++) 59 _similarities[r, c] = _similarities[c, r]; 60 _similarities[r, r] = 1; 61 for (int c = r + 1; c < _columns; c++) 62 _similarities[r, c] = (float)Kernel.KernelFunction(nodes[r], nodes[c], param); 63 } 64 } 65 66 /// <summary> 67 /// Constructor. 68 /// </summary> 69 /// <param name="rows">Nodes to use as the rows of the matrix</param> 70 /// <param name="columns">Nodes to use as the columns of the matrix</param> 71 /// <param name="param">Parameters to use when compute similarities</param> 72 public PrecomputedKernel(List<Node[]> rows, List<Node[]> columns, Parameter param) 73 { 74 _rows = rows.Count; 75 _columns = columns.Count; 76 _similarities = new float[_rows, _columns]; 77 for (int r = 0; r < _rows; r++) 78 for (int c = 0; c < _columns; c++) 79 _similarities[r, c] = (float)Kernel.KernelFunction(rows[r], columns[c], param); 44 80 } 45 81 -
trunk/sources/LibSVM/Prediction.cs
r1819 r2415 24 24 namespace SVM 25 25 { 26 /// < remarks>26 /// <summary> 27 27 /// Class containing the routines to perform class membership prediction using a trained SVM. 28 /// </ remarks>28 /// </summary> 29 29 public static class Prediction 30 30 { -
trunk/sources/LibSVM/Problem.cs
r1819 r2415 21 21 using System.IO; 22 22 using System.Collections.Generic; 23 using System.Threading; 24 using System.Globalization; 23 25 24 26 namespace SVM 25 27 { 26 /// < remarks>28 /// <summary> 27 29 /// Encapsulates a problem, or set of vectors which must be classified. 28 /// </ remarks>30 /// </summary> 29 31 [Serializable] 30 32 public class Problem … … 119 121 public static Problem Read(Stream stream) 120 122 { 123 TemporaryCulture.Start(); 124 121 125 StreamReader input = new StreamReader(stream); 122 126 List<double> vy = new List<double>(); … … 143 147 } 144 148 149 TemporaryCulture.Stop(); 150 145 151 return new Problem(vy.Count, vy.ToArray(), vx.ToArray(), max_index); 146 152 } … … 153 159 public static void Write(Stream stream, Problem problem) 154 160 { 161 TemporaryCulture.Start(); 162 155 163 StreamWriter output = new StreamWriter(stream); 156 164 for (int i = 0; i < problem.Count; i++) … … 162 170 } 163 171 output.Flush(); 172 173 TemporaryCulture.Stop(); 164 174 } 165 175 -
trunk/sources/LibSVM/Properties/AssemblyInfo.cs
r1820 r2415 6 6 // set of attributes. Change these attribute values to modify the information 7 7 // associated with an assembly. 8 [assembly: AssemblyTitle(" LibSVM")]9 [assembly: AssemblyDescription(" ")]8 [assembly: AssemblyTitle("SVM.NET")] 9 [assembly: AssemblyDescription("A .NET Support Vector Machine library adapted from libsvm")] 10 10 [assembly: AssemblyConfiguration("")] 11 11 [assembly: AssemblyCompany("Matthew Johnson")] 12 [assembly: AssemblyProduct("SVM.NET Library")]13 [assembly: AssemblyCopyright("Copyright (C) 2008 Matthew Johnson")]12 [assembly: AssemblyProduct("SVM.NET")] 13 [assembly: AssemblyCopyright("Copyright © Matthew Johnson 2009")] 14 14 [assembly: AssemblyTrademark("")] 15 15 [assembly: AssemblyCulture("")] … … 21 21 22 22 // The following GUID is for the ID of the typelib if this project is exposed to COM 23 [assembly: Guid(" 1f115aaf-e714-4341-b2ad-5e84d30b6a1a")]23 [assembly: Guid("803d4e74-70db-432e-a4da-32e7a021eec2")] 24 24 25 25 // Version information for an assembly consists of the following four values: … … 30 30 // Revision 31 31 // 32 [assembly: AssemblyVersion("2.84")] 33 [assembly: AssemblyFileVersion("2.84")] 32 // You can specify all the values or you can default the Revision and Build Numbers 33 // by using the '*' as shown below: 34 [assembly: AssemblyVersion("1.6.3")] 35 [assembly: AssemblyFileVersion("1.6.3")] -
trunk/sources/LibSVM/RangeTransform.cs
r1819 r2415 20 20 using System; 21 21 using System.IO; 22 using System.Threading; 23 using System.Globalization; 22 24 23 25 namespace SVM 24 26 { 25 /// < remarks>27 /// <summary> 26 28 /// Class which encapsulates a range transformation. 27 /// </ remarks>29 /// </summary> 28 30 public class RangeTransform : IRangeTransform 29 31 { 32 /// <summary> 33 /// Default lower bound for scaling (-1). 34 /// </summary> 35 public const int DEFAULT_LOWER_BOUND = -1; 36 /// <summary> 37 /// Default upper bound for scaling (1). 38 /// </summary> 39 public const int DEFAULT_UPPER_BOUND = 1; 40 41 /// <summary> 42 /// Determines the Range transform for the provided problem. Uses the default lower and upper bounds. 43 /// </summary> 44 /// <param name="prob">The Problem to analyze</param> 45 /// <returns>The Range transform for the problem</returns> 46 public static RangeTransform Compute(Problem prob) 47 { 48 return Compute(prob, DEFAULT_LOWER_BOUND, DEFAULT_UPPER_BOUND); 49 } 50 /// <summary> 51 /// Determines the Range transform for the provided problem. 52 /// </summary> 53 /// <param name="prob">The Problem to analyze</param> 54 /// <param name="lowerBound">The lower bound for scaling</param> 55 /// <param name="upperBound">The upper bound for scaling</param> 56 /// <returns>The Range transform for the problem</returns> 57 public static RangeTransform Compute(Problem prob, double lowerBound, double upperBound) 58 { 59 double[] minVals = new double[prob.MaxIndex]; 60 double[] maxVals = new double[prob.MaxIndex]; 61 for (int i = 0; i < prob.MaxIndex; i++) 62 { 63 minVals[i] = double.MaxValue; 64 maxVals[i] = double.MinValue; 65 } 66 for (int i = 0; i < prob.Count; i++) 67 { 68 for (int j = 0; j < prob.X[i].Length; j++) 69 { 70 int index = prob.X[i][j].Index - 1; 71 double value = prob.X[i][j].Value; 72 minVals[index] = Math.Min(minVals[index], value); 73 maxVals[index] = Math.Max(maxVals[index], value); 74 } 75 } 76 for (int i = 0; i < prob.MaxIndex; i++) 77 { 78 if (minVals[i] == double.MaxValue || maxVals[i] == double.MinValue) 79 { 80 minVals[i] = 0; 81 maxVals[i] = 0; 82 } 83 } 84 return new RangeTransform(minVals, maxVals, lowerBound, upperBound); 85 } 86 30 87 private double[] _inputStart; 31 88 private double[] _inputScale; … … 72 129 { 73 130 Node[] output = new Node[input.Length]; 74 for (int i = 0; i < _length; i++)131 for (int i = 0; i < output.Length; i++) 75 132 { 76 133 int index = input[i].Index; … … 104 161 public static void Write(Stream stream, RangeTransform r) 105 162 { 163 TemporaryCulture.Start(); 164 106 165 StreamWriter output = new StreamWriter(stream); 107 166 output.WriteLine(r._length); … … 116 175 output.WriteLine("{0} {1}", r._outputStart, r._outputScale); 117 176 output.Flush(); 177 178 TemporaryCulture.Stop(); 118 179 } 119 180 … … 161 222 public static RangeTransform Read(Stream stream) 162 223 { 224 TemporaryCulture.Start(); 225 163 226 StreamReader input = new StreamReader(stream); 164 227 int length = int.Parse(input.ReadLine()); … … 174 237 double outputStart = double.Parse(parts[0]); 175 238 double outputScale = double.Parse(parts[1]); 239 240 TemporaryCulture.Stop(); 241 176 242 return new RangeTransform(inputStart, inputScale, outputStart, outputScale, length); 177 243 } -
trunk/sources/LibSVM/Scaling.cs
r1819 r2415 22 22 namespace SVM 23 23 { 24 /// < remarks>24 /// <summary> 25 25 /// Deals with the scaling of Problems so they have uniform ranges across all dimensions in order to 26 26 /// result in better SVM performance. 27 /// </ remarks>27 /// </summary> 28 28 public static class Scaling 29 29 { 30 /// <summary>31 /// Default lower bound for scaling (-1).32 /// </summary>33 public const int DEFAULT_LOWER_BOUND = -1;34 /// <summary>35 /// Default upper bound for scaling (1).36 /// </summary>37 public const int DEFAULT_UPPER_BOUND = 1;38 39 /// <summary>40 /// Determines the Range transform for the provided problem. Uses the default lower and upper bounds.41 /// </summary>42 /// <param name="prob">The Problem to analyze</param>43 /// <returns>The Range transform for the problem</returns>44 public static RangeTransform DetermineRange(Problem prob)45 {46 return DetermineRangeTransform(prob, DEFAULT_LOWER_BOUND, DEFAULT_UPPER_BOUND);47 }48 /// <summary>49 /// Determines the Range transform for the provided problem.50 /// </summary>51 /// <param name="prob">The Problem to analyze</param>52 /// <param name="lowerBound">The lower bound for scaling</param>53 /// <param name="upperBound">The upper bound for scaling</param>54 /// <returns>The Range transform for the problem</returns>55 public static RangeTransform DetermineRangeTransform(Problem prob, double lowerBound, double upperBound)56 {57 double[] minVals = new double[prob.MaxIndex];58 double[] maxVals = new double[prob.MaxIndex];59 for (int i = 0; i < prob.MaxIndex; i++)60 {61 minVals[i] = double.MaxValue;62 maxVals[i] = double.MinValue;63 }64 for (int i = 0; i < prob.Count; i++)65 {66 for (int j = 0; j < prob.X[i].Length; j++)67 {68 int index = prob.X[i][j].Index-1;69 double value = prob.X[i][j].Value;70 minVals[index] = Math.Min(minVals[index], value);71 maxVals[index] = Math.Max(maxVals[index], value);72 }73 }74 for (int i = 0; i < prob.MaxIndex; i++)75 {76 if (minVals[i] == double.MaxValue || maxVals[i] == double.MinValue)77 {78 minVals[i] = 0;79 maxVals[i] = 0;80 }81 }82 return new RangeTransform(minVals, maxVals, lowerBound, upperBound);83 }84 30 /// <summary> 85 31 /// Scales a problem using the provided range. This will not affect the parameter. … … 88 34 /// <param name="range">The Range transform to use in scaling</param> 89 35 /// <returns>The Scaled problem</returns> 90 public static Problem Scale( Problem prob, IRangeTransform range)36 public static Problem Scale(this IRangeTransform range, Problem prob) 91 37 { 92 38 Problem scaledProblem = new Problem(prob.Count, new double[prob.Count], new Node[prob.Count][], prob.MaxIndex); -
trunk/sources/LibSVM/Solver.cs
r1819 r2415 19 19 20 20 using System; 21 using System.Linq; 21 22 using System.Collections.Generic; 22 23 using System.Diagnostics; 24 using System.IO; 23 25 24 26 namespace SVM 25 27 { 26 //27 // Kernel evaluation28 //29 // the static method k_function is for doing single kernel evaluation30 // the constructor of Kernel prepares to calculate the l*l kernel matrix31 // the member function get_Q is for getting one column from the Q Matrix32 //33 internal abstract class QMatrix34 {35 public abstract float[] get_Q(int column, int len);36 public abstract float[] get_QD();37 public abstract void swap_index(int i, int j);38 }39 40 internal abstract class Kernel : QMatrix41 {42 private Node[][] _x;43 private double[] _x_square;44 45 // Parameter46 private KernelType kernel_type;47 private int degree;48 private double gamma;49 private double coef0;50 51 public override void swap_index(int i, int j)52 {53 do { Node[] _ = _x[i]; _x[i] = _x[j]; _x[j] = _; } while (false);54 if (_x_square != null) do { double _ = _x_square[i]; _x_square[i] = _x_square[j]; _x_square[j] = _; } while (false);55 }56 57 private static double powi(double baseValue, int times)58 {59 double tmp = baseValue, ret = 1.0;60 61 for (int t = times; t > 0; t /= 2)62 {63 if (t % 2 == 1) ret *= tmp;64 tmp = tmp * tmp;65 }66 return ret;67 }68 69 private static double tanh(double x)70 {71 double e = Math.Exp(x);72 return 1.0 - 2.0 / (e * e + 1);73 }74 75 public double kernel_function(int i, int j)76 {77 switch (kernel_type)78 {79 case KernelType.LINEAR:80 return dot(_x[i], _x[j]);81 case KernelType.POLY:82 return powi(gamma * dot(_x[i], _x[j]) + coef0, degree);83 case KernelType.RBF:84 return Math.Exp(-gamma * (_x_square[i] + _x_square[j] - 2 * dot(_x[i], _x[j])));85 case KernelType.SIGMOID:86 return tanh(gamma * dot(_x[i], _x[j]) + coef0);87 case KernelType.PRECOMPUTED:88 return _x[i][(int)(_x[j][0].Value)].Value;89 default:90 return 0;91 }92 }93 94 public Kernel(int l, Node[][] x_, Parameter param)95 {96 this.kernel_type = param.KernelType;97 this.degree = param.Degree;98 this.gamma = param.Gamma;99 this.coef0 = param.Coefficient0;100 101 _x = (Node[][])x_.Clone();102 103 if (kernel_type == KernelType.RBF)104 {105 _x_square = new double[l];106 for (int i = 0; i < l; i++)107 _x_square[i] = dot(_x[i], _x[i]);108 }109 else _x_square = null;110 }111 112 public static double dot(Node[] x, Node[] y)113 {114 double sum = 0;115 int xlen = x.Length;116 int ylen = y.Length;117 int i = 0;118 int j = 0;119 while (i < xlen && j < ylen)120 {121 if (x[i].Index == y[j].Index)122 sum += x[i++].Value * y[j++].Value;123 else124 {125 if (x[i].Index > y[j].Index)126 ++j;127 else128 ++i;129 }130 }131 return sum;132 }133 134 public static double k_function(Node[] x, Node[] y, Parameter param)135 {136 switch (param.KernelType)137 {138 case KernelType.LINEAR:139 return dot(x, y);140 case KernelType.POLY:141 return powi(param.Gamma * dot(x, y) + param.Coefficient0, param.Degree);142 case KernelType.RBF:143 {144 double sum = 0;145 int xlen = x.Length;146 int ylen = y.Length;147 int i = 0;148 int j = 0;149 while (i < xlen && j < ylen)150 {151 if (x[i].Index == y[j].Index)152 {153 double d = x[i++].Value - y[j++].Value;154 sum += d * d;155 }156 else if (x[i].Index > y[j].Index)157 {158 sum += y[j].Value * y[j].Value;159 ++j;160 }161 else162 {163 sum += x[i].Value * x[i].Value;164 ++i;165 }166 }167 168 while (i < xlen)169 {170 sum += x[i].Value * x[i].Value;171 ++i;172 }173 174 while (j < ylen)175 {176 sum += y[j].Value * y[j].Value;177 ++j;178 }179 180 return Math.Exp(-param.Gamma * sum);181 }182 case KernelType.SIGMOID:183 return tanh(param.Gamma * dot(x, y) + param.Coefficient0);184 case KernelType.PRECOMPUTED:185 return x[(int)(y[0].Value)].Value;186 default:187 return 0;188 }189 }190 }191 28 192 29 // An SMO algorithm in Fan et al., JMLR 6(2005), p. 1889--1918 193 30 // Solves: 194 31 // 195 // min 0.5(\alpha^T Q \alpha) + p^T \alpha32 // Min 0.5(\alpha^T Q \alpha) + p^T \alpha 196 33 // 197 34 // y^T \alpha = \delta … … 211 48 { 212 49 protected int active_size; 213 protected s hort[] y;50 protected sbyte[] y; 214 51 protected double[] G; // gradient of objective function 215 pr otectedconst byte LOWER_BOUND = 0;216 pr otectedconst byte UPPER_BOUND = 1;217 pr otectedconst byte FREE = 2;218 pr otectedbyte[] alpha_status; // LOWER_BOUND, UPPER_BOUND, FREE219 pr otecteddouble[] alpha;220 protected QMatrix Q;52 private const byte LOWER_BOUND = 0; 53 private const byte UPPER_BOUND = 1; 54 private const byte FREE = 2; 55 private byte[] alpha_status; // LOWER_BOUND, UPPER_BOUND, FREE 56 private double[] alpha; 57 protected IQMatrix Q; 221 58 protected float[] QD; 222 protected double eps;223 pr otecteddouble Cp, Cn;224 pr otecteddouble[] p;225 pr otectedint[] active_set;226 pr otecteddouble[] G_bar; // gradient, if we treat free variables as 059 protected double EPS; 60 private double Cp, Cn; 61 private double[] p; 62 private int[] active_set; 63 private double[] G_bar; // gradient, if we treat free variables as 0 227 64 protected int l; 228 protected bool unshrink ed; // XXX65 protected bool unshrink; // XXX 229 66 230 67 protected const double INF = double.PositiveInfinity; 231 68 232 pr otecteddouble get_C(int i)69 private double get_C(int i) 233 70 { 234 71 return (y[i] > 0) ? Cp : Cn; 235 72 } 236 protected void update_alpha_status(int i) 73 74 private void update_alpha_status(int i) 237 75 { 238 76 if (alpha[i] >= get_C(i)) … … 242 80 else alpha_status[i] = FREE; 243 81 } 82 244 83 protected bool is_upper_bound(int i) { return alpha_status[i] == UPPER_BOUND; } 245 84 protected bool is_lower_bound(int i) { return alpha_status[i] == LOWER_BOUND; } 246 protected bool is_free(int i) { return alpha_status[i] == FREE; } 247 248 // java: information about solution except alpha, 249 // because we cannot return multiple values otherwise... 250 internal class SolutionInfo 85 86 private bool is_free(int i) { return alpha_status[i] == FREE; } 87 88 public class SolutionInfo 251 89 { 252 90 public double obj; … … 259 97 protected void swap_index(int i, int j) 260 98 { 261 Q. swap_index(i, j);262 do { short _ = y[i]; y[i] = y[j]; y[j] = _; } while (false);263 do { double _ = G[i]; G[i] = G[j]; G[j] = _; } while (false);264 do { byte _ = alpha_status[i]; alpha_status[i] = alpha_status[j]; alpha_status[j] = _; } while (false);265 do { double _ = alpha[i]; alpha[i] = alpha[j]; alpha[j] = _; } while (false);266 do { double _ = p[i]; p[i] = p[j]; p[j] = _; } while (false);267 do { int _ = active_set[i]; active_set[i] = active_set[j]; active_set[j] = _; } while (false);268 do { double _ = G_bar[i]; G_bar[i] = G_bar[j]; G_bar[j] = _; } while (false);99 Q.SwapIndex(i, j); 100 y.SwapIndex(i, j); 101 G.SwapIndex(i, j); 102 alpha_status.SwapIndex(i, j); 103 alpha.SwapIndex(i, j); 104 p.SwapIndex(i, j); 105 active_set.SwapIndex(i, j); 106 G_bar.SwapIndex(i, j); 269 107 } 270 108 … … 275 113 if (active_size == l) return; 276 114 277 int i; 278 for (i = active_size; i < l; i++) 279 G[i] = G_bar[i] + p[i]; 280 281 for (i = 0; i < active_size; i++) 282 if (is_free(i)) 283 { 284 float[] Q_i = Q.get_Q(i, l); 285 double alpha_i = alpha[i]; 286 for (int j = active_size; j < l; j++) 287 G[j] += alpha_i * Q_i[j]; 288 } 289 } 290 291 public virtual void Solve(int l, QMatrix Q, double[] p_, short[] y_, 292 double[] alpha_, double Cp, double Cn, double eps, SolutionInfo si, bool shrinking) 115 int i, j; 116 int nr_free = 0; 117 118 for (j = active_size; j < l; j++) 119 G[j] = G_bar[j] + p[j]; 120 121 for (j = 0; j < active_size; j++) 122 if (is_free(j)) 123 nr_free++; 124 125 if (2 * nr_free < active_size) 126 Procedures.info("\nWarning: using -h 0 may be faster\n"); 127 128 if (nr_free * l > 2 * active_size * (l - active_size)) 129 { 130 for (i = active_size; i < l; i++) 131 { 132 float[] Q_i = Q.GetQ(i, active_size); 133 for (j = 0; j < active_size; j++) 134 if (is_free(j)) 135 G[i] += alpha[j] * Q_i[j]; 136 } 137 } 138 else 139 { 140 for (i = 0; i < active_size; i++) 141 if (is_free(i)) 142 { 143 float[] Q_i = Q.GetQ(i, l); 144 double alpha_i = alpha[i]; 145 for (j = active_size; j < l; j++) 146 G[j] += alpha_i * Q_i[j]; 147 } 148 } 149 } 150 151 public virtual void Solve(int l, IQMatrix Q, double[] p_, sbyte[] y_, double[] alpha_, double Cp, double Cn, double eps, SolutionInfo si, bool shrinking) 293 152 { 294 153 this.l = l; 295 154 this.Q = Q; 296 QD = Q. get_QD();155 QD = Q.GetQD(); 297 156 p = (double[])p_.Clone(); 298 y = (s hort[])y_.Clone();157 y = (sbyte[])y_.Clone(); 299 158 alpha = (double[])alpha_.Clone(); 300 159 this.Cp = Cp; 301 160 this.Cn = Cn; 302 this. eps= eps;303 this.unshrink ed= false;161 this.EPS = eps; 162 this.unshrink = false; 304 163 305 164 // initialize alpha_status … … 331 190 if (!is_lower_bound(i)) 332 191 { 333 float[] Q_i = Q. get_Q(i, l);192 float[] Q_i = Q.GetQ(i, l); 334 193 double alpha_i = alpha[i]; 335 194 int j; … … 356 215 counter = Math.Min(l, 1000); 357 216 if (shrinking) do_shrinking(); 358 Debug.Write(".");217 Procedures.info("."); 359 218 } 360 219 … … 365 224 // reset active set size and check 366 225 active_size = l; 367 Debug.Write("*");226 Procedures.info("*"); 368 227 if (select_working_set(working_set) != 0) 369 228 break; … … 379 238 // update alpha[i] and alpha[j], handle bounds carefully 380 239 381 float[] Q_i = Q. get_Q(i, active_size);382 float[] Q_j = Q. get_Q(j, active_size);240 float[] Q_i = Q.GetQ(i, active_size); 241 float[] Q_j = Q.GetQ(j, active_size); 383 242 384 243 double C_i = get_C(i); … … 495 354 if (ui != is_upper_bound(i)) 496 355 { 497 Q_i = Q. get_Q(i, l);356 Q_i = Q.GetQ(i, l); 498 357 if (ui) 499 358 for (k = 0; k < l; k++) … … 506 365 if (uj != is_upper_bound(j)) 507 366 { 508 Q_j = Q. get_Q(j, l);367 Q_j = Q.GetQ(j, l); 509 368 if (uj) 510 369 for (k = 0; k < l; k++) … … 541 400 si.upper_bound_n = Cn; 542 401 543 Debug.Write("\noptimization finished, #iter = " + iter + "\n");402 Procedures.info("\noptimization finished, #iter = " + iter + "\n"); 544 403 } 545 404 546 405 // return 1 if already optimal, return 0 otherwise 547 protected virtualint select_working_set(int[] working_set)406 int select_working_set(int[] working_set) 548 407 { 549 408 // return i,j such that 550 // i: maximizes -y_i * grad(f)_i, i in I_up(\alpha)409 // i: Maximizes -y_i * grad(f)_i, i in I_up(\alpha) 551 410 // j: mimimizes the decrease of obj value 552 411 // (if quadratic coefficeint <= 0, replace it with tau) 553 412 // -y_j*grad(f)_j < -y_i*grad(f)_i, j in I_low(\alpha) 554 413 555 double G max = -INF;556 double G max2 = -INF;557 int G max_idx = -1;558 int G min_idx = -1;559 double obj_diff_ min = INF;414 double GMax = -INF; 415 double GMax2 = -INF; 416 int GMax_idx = -1; 417 int GMin_idx = -1; 418 double obj_diff_Min = INF; 560 419 561 420 for (int t = 0; t < active_size; t++) … … 563 422 { 564 423 if (!is_upper_bound(t)) 565 if (-G[t] >= G max)566 { 567 G max = -G[t];568 G max_idx = t;424 if (-G[t] >= GMax) 425 { 426 GMax = -G[t]; 427 GMax_idx = t; 569 428 } 570 429 } … … 572 431 { 573 432 if (!is_lower_bound(t)) 574 if (G[t] >= G max)575 { 576 G max = G[t];577 G max_idx = t;578 } 579 } 580 581 int i = G max_idx;433 if (G[t] >= GMax) 434 { 435 GMax = G[t]; 436 GMax_idx = t; 437 } 438 } 439 440 int i = GMax_idx; 582 441 float[] Q_i = null; 583 if (i != -1) // null Q_i not accessed: G max=-INF if i=-1584 Q_i = Q. get_Q(i, active_size);442 if (i != -1) // null Q_i not accessed: GMax=-INF if i=-1 443 Q_i = Q.GetQ(i, active_size); 585 444 586 445 for (int j = 0; j < active_size; j++) … … 590 449 if (!is_lower_bound(j)) 591 450 { 592 double grad_diff = G max + G[j];593 if (G[j] >= G max2)594 G max2 = G[j];451 double grad_diff = GMax + G[j]; 452 if (G[j] >= GMax2) 453 GMax2 = G[j]; 595 454 if (grad_diff > 0) 596 455 { 597 456 double obj_diff; 598 double quad_coef = Q_i[i] + QD[j] - 2 * y[i] * Q_i[j];457 double quad_coef = Q_i[i] + QD[j] - 2.0 * y[i] * Q_i[j]; 599 458 if (quad_coef > 0) 600 459 obj_diff = -(grad_diff * grad_diff) / quad_coef; … … 602 461 obj_diff = -(grad_diff * grad_diff) / 1e-12; 603 462 604 if (obj_diff <= obj_diff_ min)463 if (obj_diff <= obj_diff_Min) 605 464 { 606 G min_idx = j;607 obj_diff_ min = obj_diff;465 GMin_idx = j; 466 obj_diff_Min = obj_diff; 608 467 } 609 468 } … … 614 473 if (!is_upper_bound(j)) 615 474 { 616 double grad_diff = G max - G[j];617 if (-G[j] >= G max2)618 G max2 = -G[j];475 double grad_diff = GMax - G[j]; 476 if (-G[j] >= GMax2) 477 GMax2 = -G[j]; 619 478 if (grad_diff > 0) 620 479 { 621 480 double obj_diff; 622 double quad_coef = Q_i[i] + QD[j] + 2 * y[i] * Q_i[j];481 double quad_coef = Q_i[i] + QD[j] + 2.0 * y[i] * Q_i[j]; 623 482 if (quad_coef > 0) 624 483 obj_diff = -(grad_diff * grad_diff) / quad_coef; … … 626 485 obj_diff = -(grad_diff * grad_diff) / 1e-12; 627 486 628 if (obj_diff <= obj_diff_ min)487 if (obj_diff <= obj_diff_Min) 629 488 { 630 G min_idx = j;631 obj_diff_ min = obj_diff;489 GMin_idx = j; 490 obj_diff_Min = obj_diff; 632 491 } 633 492 } … … 636 495 } 637 496 638 if (G max + Gmax2 < eps)497 if (GMax + GMax2 < EPS) 639 498 return 1; 640 499 641 working_set[0] = G max_idx;642 working_set[1] = G min_idx;500 working_set[0] = GMax_idx; 501 working_set[1] = GMin_idx; 643 502 return 0; 644 503 } 645 504 646 private bool be_shrunk en(int i, double Gmax1, double Gmax2)505 private bool be_shrunk(int i, double GMax1, double GMax2) 647 506 { 648 507 if (is_upper_bound(i)) 649 508 { 650 509 if (y[i] == +1) 651 return (-G[i] > G max1);510 return (-G[i] > GMax1); 652 511 else 653 return (-G[i] > G max2);512 return (-G[i] > GMax2); 654 513 } 655 514 else if (is_lower_bound(i)) 656 515 { 657 516 if (y[i] == +1) 658 return (G[i] > G max2);517 return (G[i] > GMax2); 659 518 else 660 return (G[i] > G max1);519 return (G[i] > GMax1); 661 520 } 662 521 else … … 664 523 } 665 524 666 protected virtualvoid do_shrinking()525 void do_shrinking() 667 526 { 668 527 int i; 669 double G max1 = -INF; // max { -y_i * grad(f)_i | i in I_up(\alpha) }670 double G max2 = -INF; // max { y_i * grad(f)_i | i in I_low(\alpha) }671 672 // find maximal violating pair first528 double GMax1 = -INF; // Max { -y_i * grad(f)_i | i in I_up(\alpha) } 529 double GMax2 = -INF; // Max { y_i * grad(f)_i | i in I_low(\alpha) } 530 531 // find Maximal violating pair first 673 532 for (i = 0; i < active_size; i++) 674 533 { … … 677 536 if (!is_upper_bound(i)) 678 537 { 679 if (-G[i] >= G max1)680 G max1 = -G[i];538 if (-G[i] >= GMax1) 539 GMax1 = -G[i]; 681 540 } 682 541 if (!is_lower_bound(i)) 683 542 { 684 if (G[i] >= G max2)685 G max2 = G[i];543 if (G[i] >= GMax2) 544 GMax2 = G[i]; 686 545 } 687 546 } … … 690 549 if (!is_upper_bound(i)) 691 550 { 692 if (-G[i] >= G max2)693 G max2 = -G[i];551 if (-G[i] >= GMax2) 552 GMax2 = -G[i]; 694 553 } 695 554 if (!is_lower_bound(i)) 696 555 { 697 if (G[i] >= Gmax1) 698 Gmax1 = G[i]; 699 } 700 } 701 } 702 703 // shrink 556 if (G[i] >= GMax1) 557 GMax1 = G[i]; 558 } 559 } 560 } 561 562 if (unshrink == false && GMax1 + GMax2 <= EPS * 10) 563 { 564 unshrink = true; 565 reconstruct_gradient(); 566 active_size = l; 567 } 704 568 705 569 for (i = 0; i < active_size; i++) 706 if (be_shrunk en(i, Gmax1, Gmax2))570 if (be_shrunk(i, GMax1, GMax2)) 707 571 { 708 572 active_size--; 709 573 while (active_size > i) 710 574 { 711 if (!be_shrunk en(active_size, Gmax1, Gmax2))575 if (!be_shrunk(active_size, GMax1, GMax2)) 712 576 { 713 577 swap_index(i, active_size); … … 717 581 } 718 582 } 719 720 // unshrink, check all variables again before sealed iterations 721 722 if (unshrinked || Gmax1 + Gmax2 > eps * 10) return; 723 724 unshrinked = true; 725 reconstruct_gradient(); 726 727 for (i = l - 1; i >= active_size; i--) 728 if (!be_shrunken(i, Gmax1, Gmax2)) 729 { 730 while (active_size < i) 731 { 732 if (be_shrunken(active_size, Gmax1, Gmax2)) 733 { 734 swap_index(i, active_size); 735 break; 736 } 737 active_size++; 738 } 739 active_size++; 740 } 741 } 742 743 protected virtual double calculate_rho() 583 } 584 585 double calculate_rho() 744 586 { 745 587 double r; … … 786 628 // additional constraint: e^T \alpha = constant 787 629 // 788 sealedclass Solver_NU : Solver630 class Solver_NU : Solver 789 631 { 790 632 private SolutionInfo si; 791 633 792 public override void Solve(int l, QMatrix Q, double[] p, short[] y,634 public sealed override void Solve(int l, IQMatrix Q, double[] p, sbyte[] y, 793 635 double[] alpha, double Cp, double Cn, double eps, 794 636 SolutionInfo si, bool shrinking) … … 799 641 800 642 // return 1 if already optimal, return 0 otherwise 801 pr otected override int select_working_set(int[] working_set)643 private int select_working_set(int[] working_set) 802 644 { 803 645 // return i,j such that y_i = y_j and 804 // i: maximizes -y_i * grad(f)_i, i in I_up(\alpha)805 // j: minimizes the decrease of obj value646 // i: Maximizes -y_i * grad(f)_i, i in I_up(\alpha) 647 // j: Minimizes the decrease of obj value 806 648 // (if quadratic coefficeint <= 0, replace it with tau) 807 649 // -y_j*grad(f)_j < -y_i*grad(f)_i, j in I_low(\alpha) 808 650 809 double G maxp = -INF;810 double G maxp2 = -INF;811 int G maxp_idx = -1;812 813 double G maxn = -INF;814 double G maxn2 = -INF;815 int G maxn_idx = -1;816 817 int G min_idx = -1;818 double obj_diff_ min = INF;651 double GMaxp = -INF; 652 double GMaxp2 = -INF; 653 int GMaxp_idx = -1; 654 655 double GMaxn = -INF; 656 double GMaxn2 = -INF; 657 int GMaxn_idx = -1; 658 659 int GMin_idx = -1; 660 double obj_diff_Min = INF; 819 661 820 662 for (int t = 0; t < active_size; t++) … … 822 664 { 823 665 if (!is_upper_bound(t)) 824 if (-G[t] >= G maxp)825 { 826 G maxp = -G[t];827 G maxp_idx = t;666 if (-G[t] >= GMaxp) 667 { 668 GMaxp = -G[t]; 669 GMaxp_idx = t; 828 670 } 829 671 } … … 831 673 { 832 674 if (!is_lower_bound(t)) 833 if (G[t] >= G maxn)834 { 835 G maxn = G[t];836 G maxn_idx = t;837 } 838 } 839 840 int ip = G maxp_idx;841 int iN = G maxn_idx;675 if (G[t] >= GMaxn) 676 { 677 GMaxn = G[t]; 678 GMaxn_idx = t; 679 } 680 } 681 682 int ip = GMaxp_idx; 683 int iN = GMaxn_idx; 842 684 float[] Q_ip = null; 843 685 float[] Q_in = null; 844 if (ip != -1) // null Q_ip not accessed: G maxp=-INF if ip=-1845 Q_ip = Q. get_Q(ip, active_size);686 if (ip != -1) // null Q_ip not accessed: GMaxp=-INF if ip=-1 687 Q_ip = Q.GetQ(ip, active_size); 846 688 if (iN != -1) 847 Q_in = Q. get_Q(iN, active_size);689 Q_in = Q.GetQ(iN, active_size); 848 690 849 691 for (int j = 0; j < active_size; j++) … … 853 695 if (!is_lower_bound(j)) 854 696 { 855 double grad_diff = G maxp + G[j];856 if (G[j] >= G maxp2)857 G maxp2 = G[j];697 double grad_diff = GMaxp + G[j]; 698 if (G[j] >= GMaxp2) 699 GMaxp2 = G[j]; 858 700 if (grad_diff > 0) 859 701 { … … 865 707 obj_diff = -(grad_diff * grad_diff) / 1e-12; 866 708 867 if (obj_diff <= obj_diff_ min)709 if (obj_diff <= obj_diff_Min) 868 710 { 869 G min_idx = j;870 obj_diff_ min = obj_diff;711 GMin_idx = j; 712 obj_diff_Min = obj_diff; 871 713 } 872 714 } … … 877 719 if (!is_upper_bound(j)) 878 720 { 879 double grad_diff = G maxn - G[j];880 if (-G[j] >= G maxn2)881 G maxn2 = -G[j];721 double grad_diff = GMaxn - G[j]; 722 if (-G[j] >= GMaxn2) 723 GMaxn2 = -G[j]; 882 724 if (grad_diff > 0) 883 725 { … … 889 731 obj_diff = -(grad_diff * grad_diff) / 1e-12; 890 732 891 if (obj_diff <= obj_diff_ min)733 if (obj_diff <= obj_diff_Min) 892 734 { 893 G min_idx = j;894 obj_diff_ min = obj_diff;735 GMin_idx = j; 736 obj_diff_Min = obj_diff; 895 737 } 896 738 } … … 899 741 } 900 742 901 if (Math.Max(G maxp + Gmaxp2, Gmaxn + Gmaxn2) < eps)743 if (Math.Max(GMaxp + GMaxp2, GMaxn + GMaxn2) < EPS) 902 744 return 1; 903 745 904 if (y[G min_idx] == +1)905 working_set[0] = G maxp_idx;746 if (y[GMin_idx] == +1) 747 working_set[0] = GMaxp_idx; 906 748 else 907 working_set[0] = G maxn_idx;908 working_set[1] = G min_idx;749 working_set[0] = GMaxn_idx; 750 working_set[1] = GMin_idx; 909 751 910 752 return 0; 911 753 } 912 754 913 private bool be_shrunk en(int i, double Gmax1, double Gmax2, double Gmax3, double Gmax4)755 private bool be_shrunk(int i, double GMax1, double GMax2, double GMax3, double GMax4) 914 756 { 915 757 if (is_upper_bound(i)) 916 758 { 917 759 if (y[i] == +1) 918 return (-G[i] > G max1);760 return (-G[i] > GMax1); 919 761 else 920 return (-G[i] > G max4);762 return (-G[i] > GMax4); 921 763 } 922 764 else if (is_lower_bound(i)) 923 765 { 924 766 if (y[i] == +1) 925 return (G[i] > G max2);767 return (G[i] > GMax2); 926 768 else 927 return (G[i] > G max3);769 return (G[i] > GMax3); 928 770 } 929 771 else … … 931 773 } 932 774 933 pr otected override void do_shrinking()934 { 935 double G max1 = -INF; // max { -y_i * grad(f)_i | y_i = +1, i in I_up(\alpha) }936 double G max2 = -INF; // max { y_i * grad(f)_i | y_i = +1, i in I_low(\alpha) }937 double G max3 = -INF; // max { -y_i * grad(f)_i | y_i = -1, i in I_up(\alpha) }938 double G max4 = -INF; // max { y_i * grad(f)_i | y_i = -1, i in I_low(\alpha) }939 940 // find maximal violating pair first775 private void do_shrinking() 776 { 777 double GMax1 = -INF; // Max { -y_i * grad(f)_i | y_i = +1, i in I_up(\alpha) } 778 double GMax2 = -INF; // Max { y_i * grad(f)_i | y_i = +1, i in I_low(\alpha) } 779 double GMax3 = -INF; // Max { -y_i * grad(f)_i | y_i = -1, i in I_up(\alpha) } 780 double GMax4 = -INF; // Max { y_i * grad(f)_i | y_i = -1, i in I_low(\alpha) } 781 782 // find Maximal violating pair first 941 783 int i; 942 784 for (i = 0; i < active_size; i++) … … 946 788 if (y[i] == +1) 947 789 { 948 if (-G[i] > G max1) Gmax1 = -G[i];949 } 950 else if (-G[i] > G max4) Gmax4 = -G[i];790 if (-G[i] > GMax1) GMax1 = -G[i]; 791 } 792 else if (-G[i] > GMax4) GMax4 = -G[i]; 951 793 } 952 794 if (!is_lower_bound(i)) … … 954 796 if (y[i] == +1) 955 797 { 956 if (G[i] > Gmax2) Gmax2 = G[i]; 957 } 958 else if (G[i] > Gmax3) Gmax3 = G[i]; 959 } 960 } 961 962 // shrinking 798 if (G[i] > GMax2) GMax2 = G[i]; 799 } 800 else if (G[i] > GMax3) GMax3 = G[i]; 801 } 802 } 803 804 if (unshrink == false && Math.Max(GMax1 + GMax2, GMax3 + GMax4) <= EPS * 10) 805 { 806 unshrink = true; 807 reconstruct_gradient(); 808 active_size = l; 809 } 963 810 964 811 for (i = 0; i < active_size; i++) 965 if (be_shrunk en(i, Gmax1, Gmax2, Gmax3, Gmax4))812 if (be_shrunk(i, GMax1, GMax2, GMax3, GMax4)) 966 813 { 967 814 active_size--; 968 815 while (active_size > i) 969 816 { 970 if (!be_shrunk en(active_size, Gmax1, Gmax2, Gmax3, Gmax4))817 if (!be_shrunk(active_size, GMax1, GMax2, GMax3, GMax4)) 971 818 { 972 819 swap_index(i, active_size); … … 976 823 } 977 824 } 978 979 if (unshrinked || Math.Max(Gmax1 + Gmax2, Gmax3 + Gmax4) > eps * 10) return; 980 981 unshrinked = true; 982 reconstruct_gradient(); 983 984 for (i = l - 1; i >= active_size; i--) 985 if (!be_shrunken(i, Gmax1, Gmax2, Gmax3, Gmax4)) 986 { 987 while (active_size < i) 988 { 989 if (be_shrunken(active_size, Gmax1, Gmax2, Gmax3, Gmax4)) 990 { 991 swap_index(i, active_size); 992 break; 993 } 994 active_size++; 995 } 996 active_size++; 997 } 998 } 999 1000 protected override double calculate_rho() 825 } 826 827 private double calculate_rho() 1001 828 { 1002 829 int nr_free1 = 0, nr_free2 = 0; … … 1054 881 class SVC_Q : Kernel 1055 882 { 1056 private s hort[] y;883 private sbyte[] y; 1057 884 private Cache cache; 1058 885 private float[] QD; 1059 886 1060 public SVC_Q(Problem prob, Parameter param, s hort[] y_) : base(prob.Count, prob.X, param)1061 { 1062 y = (s hort[])y_.Clone();887 public SVC_Q(Problem prob, Parameter param, sbyte[] y_) : base(prob.Count, prob.X, param) 888 { 889 y = (sbyte[])y_.Clone(); 1063 890 cache = new Cache(prob.Count, (long)(param.CacheSize * (1 << 20))); 1064 891 QD = new float[prob.Count]; 1065 892 for (int i = 0; i < prob.Count; i++) 1066 QD[i] = (float) kernel_function(i, i);1067 } 1068 1069 public override float[] get_Q(int i, int len)1070 { 1071 float[] [] data = new float[1][];1072 int start ;1073 if ((start = cache. get_data(i,data, len)) < len)1074 { 1075 for ( intj = start; j < len; j++)1076 data[ 0][j] = (float)(y[i] * y[j] * kernel_function(i, j));1077 } 1078 return data [0];1079 } 1080 1081 public override float[] get_QD()893 QD[i] = (float)KernelFunction(i, i); 894 } 895 896 public override sealed float[] GetQ(int i, int len) 897 { 898 float[] data = null; 899 int start, j; 900 if ((start = cache.GetData(i, ref data, len)) < len) 901 { 902 for (j = start; j < len; j++) 903 data[j] = (float)(y[i] * y[j] * KernelFunction(i, j)); 904 } 905 return data; 906 } 907 908 public override sealed float[] GetQD() 1082 909 { 1083 910 return QD; 1084 911 } 1085 912 1086 public override void swap_index(int i, int j)1087 { 1088 cache. swap_index(i, j);1089 base. swap_index(i, j);1090 do { short _ = y[i]; y[i] = y[j]; y[j] = _; } while (false);1091 do { float _ = QD[i]; QD[i] = QD[j]; QD[j] = _; } while (false);913 public override sealed void SwapIndex(int i, int j) 914 { 915 cache.SwapIndex(i, j); 916 base.SwapIndex(i, j); 917 y.SwapIndex(i, j); 918 QD.SwapIndex(i, j); 1092 919 } 1093 920 } … … 1098 925 private float[] QD; 1099 926 1100 public ONE_CLASS_Q(Problem prob, Parameter param) : base(prob.Count, prob.X, param)927 public ONE_CLASS_Q(Problem prob, Parameter param) : base(prob.Count, prob.X, param) 1101 928 { 1102 929 cache = new Cache(prob.Count, (long)(param.CacheSize * (1 << 20))); 1103 930 QD = new float[prob.Count]; 1104 931 for (int i = 0; i < prob.Count; i++) 1105 QD[i] = (float) kernel_function(i, i);1106 } 1107 1108 public override float[] get_Q(int i, int len)1109 { 1110 float[] [] data = new float[1][];1111 int start ;1112 if ((start = cache. get_data(i,data, len)) < len)1113 { 1114 for ( intj = start; j < len; j++)1115 data[ 0][j] = (float)kernel_function(i, j);1116 } 1117 return data [0];1118 } 1119 1120 public override float[] get_QD()932 QD[i] = (float)KernelFunction(i, i); 933 } 934 935 public override sealed float[] GetQ(int i, int len) 936 { 937 float[] data = null; 938 int start, j; 939 if ((start = cache.GetData(i, ref data, len)) < len) 940 { 941 for (j = start; j < len; j++) 942 data[j] = (float)KernelFunction(i, j); 943 } 944 return data; 945 } 946 947 public override sealed float[] GetQD() 1121 948 { 1122 949 return QD; 1123 950 } 1124 951 1125 public override void swap_index(int i, int j)1126 { 1127 cache. swap_index(i, j);1128 base. swap_index(i, j);1129 do { float _ = QD[i]; QD[i] = QD[j]; QD[j] = _; } while (false);952 public override sealed void SwapIndex(int i, int j) 953 { 954 cache.SwapIndex(i, j); 955 base.SwapIndex(i, j); 956 QD.SwapIndex(i, j); 1130 957 } 1131 958 } … … 1135 962 private int l; 1136 963 private Cache cache; 1137 private s hort[] sign;964 private sbyte[] sign; 1138 965 private int[] index; 1139 966 private int next_buffer; … … 1141 968 private float[] QD; 1142 969 1143 public SVR_Q(Problem prob, Parameter param) 1144 : base(prob.Count, prob.X, param) 970 public SVR_Q(Problem prob, Parameter param) : base(prob.Count, prob.X, param) 1145 971 { 1146 972 l = prob.Count; 1147 973 cache = new Cache(l, (long)(param.CacheSize * (1 << 20))); 1148 974 QD = new float[2 * l]; 1149 sign = new s hort[2 * l];975 sign = new sbyte[2 * l]; 1150 976 index = new int[2 * l]; 1151 977 for (int k = 0; k < l; k++) … … 1155 981 index[k] = k; 1156 982 index[k + l] = k; 1157 QD[k] = (float) kernel_function(k, k);983 QD[k] = (float)KernelFunction(k, k); 1158 984 QD[k + l] = QD[k]; 1159 985 } … … 1164 990 } 1165 991 1166 public override void swap_index(int i, int j)1167 { 1168 do { short _ = sign[i]; sign[i] = sign[j]; sign[j] = _; } while (false);1169 do { int _ = index[i]; index[i] = index[j]; index[j] = _; } while (false);1170 do { float _ = QD[i]; QD[i] = QD[j]; QD[j] = _; } while (false);1171 } 1172 1173 public override float[] get_Q(int i, int len)1174 { 1175 float[] [] data = new float[1][];1176 int real_i = index[i];1177 if (cache. get_data(real_i,data, l) < l)1178 { 1179 for ( intj = 0; j < l; j++)1180 data[ 0][j] = (float)kernel_function(real_i, j);992 public override sealed void SwapIndex(int i, int j) 993 { 994 sign.SwapIndex(i, j); 995 index.SwapIndex(i, j); 996 QD.SwapIndex(i, j); 997 } 998 999 public override sealed float[] GetQ(int i, int len) 1000 { 1001 float[] data = null; 1002 int j, real_i = index[i]; 1003 if (cache.GetData(real_i, ref data, l) < l) 1004 { 1005 for (j = 0; j < l; j++) 1006 data[j] = (float)KernelFunction(real_i, j); 1181 1007 } 1182 1008 … … 1184 1010 float[] buf = buffer[next_buffer]; 1185 1011 next_buffer = 1 - next_buffer; 1186 s hortsi = sign[i];1187 for ( intj = 0; j < len; j++)1188 buf[j] = si * sign[j] * data[0][index[j]];1012 sbyte si = sign[i]; 1013 for (j = 0; j < len; j++) 1014 buf[j] = (float)si * sign[j] * data[index[j]]; 1189 1015 return buf; 1190 1016 } 1191 1017 1192 public override float[] get_QD()1018 public override sealed float[] GetQD() 1193 1019 { 1194 1020 return QD; … … 1196 1022 } 1197 1023 1198 internal staticclass Procedures1024 internal class Procedures 1199 1025 { 1026 private static bool _verbose; 1027 public static bool IsVerbose 1028 { 1029 get 1030 { 1031 return _verbose; 1032 } 1033 set 1034 { 1035 _verbose = value; 1036 } 1037 } 1200 1038 // 1201 1039 // construct and solve various formulations 1202 1040 // 1041 public const int LIBSVM_VERSION = 289; 1042 1043 public static TextWriter svm_print_string = Console.Out; 1044 1045 public static void info(string s) 1046 { 1047 if(_verbose) 1048 svm_print_string.Write(s); 1049 } 1050 1203 1051 private static void solve_c_svc(Problem prob, Parameter param, 1204 1052 double[] alpha, Solver.SolutionInfo si, … … 1206 1054 { 1207 1055 int l = prob.Count; 1208 double[] minus_ones = new double[l];1209 s hort[] y = new short[l];1056 double[] Minus_ones = new double[l]; 1057 sbyte[] y = new sbyte[l]; 1210 1058 1211 1059 int i; … … 1214 1062 { 1215 1063 alpha[i] = 0; 1216 minus_ones[i] = -1;1064 Minus_ones[i] = -1; 1217 1065 if (prob.Y[i] > 0) y[i] = +1; else y[i] = -1; 1218 1066 } 1219 1067 1220 1068 Solver s = new Solver(); 1221 s.Solve(l, new SVC_Q(prob, param, y), minus_ones, y,1069 s.Solve(l, new SVC_Q(prob, param, y), Minus_ones, y, 1222 1070 alpha, Cp, Cn, param.EPS, si, param.Shrinking); 1223 1071 … … 1227 1075 1228 1076 if (Cp == Cn) 1229 Debug.Write("nu = " + sum_alpha / (Cp * prob.Count) + "\n");1077 Procedures.info("nu = " + sum_alpha / (Cp * prob.Count) + "\n"); 1230 1078 1231 1079 for (i = 0; i < l; i++) … … 1240 1088 double nu = param.Nu; 1241 1089 1242 s hort[] y = new short[l];1090 sbyte[] y = new sbyte[l]; 1243 1091 1244 1092 for (i = 0; i < l; i++) … … 1269 1117 1270 1118 Solver_NU s = new Solver_NU(); 1271 s.Solve(l, new SVC_Q(prob, param, y), zeros, y, 1272 alpha, 1.0, 1.0, param.EPS, si, param.Shrinking); 1119 s.Solve(l, new SVC_Q(prob, param, y), zeros, y, alpha, 1.0, 1.0, param.EPS, si, param.Shrinking); 1273 1120 double r = si.r; 1274 1121 1275 Debug.Write("C = " + 1 / r + "\n");1122 Procedures.info("C = " + 1 / r + "\n"); 1276 1123 1277 1124 for (i = 0; i < l; i++) … … 1285 1132 1286 1133 private static void solve_one_class(Problem prob, Parameter param, 1287 1134 double[] alpha, Solver.SolutionInfo si) 1288 1135 { 1289 1136 int l = prob.Count; 1290 1137 double[] zeros = new double[l]; 1291 s hort[] ones = new short[l];1138 sbyte[] ones = new sbyte[l]; 1292 1139 int i; 1293 1140 … … 1308 1155 1309 1156 Solver s = new Solver(); 1310 s.Solve(l, new ONE_CLASS_Q(prob, param), zeros, ones, 1311 alpha, 1.0, 1.0, param.EPS, si, param.Shrinking); 1312 } 1313 1314 private static void solve_epsilon_svr(Problem prob, Parameter param, 1315 double[] alpha, Solver.SolutionInfo si) 1157 s.Solve(l, new ONE_CLASS_Q(prob, param), zeros, ones, alpha, 1.0, 1.0, param.EPS, si, param.Shrinking); 1158 } 1159 1160 private static void solve_epsilon_svr(Problem prob, Parameter param, double[] alpha, Solver.SolutionInfo si) 1316 1161 { 1317 1162 int l = prob.Count; 1318 1163 double[] alpha2 = new double[2 * l]; 1319 1164 double[] linear_term = new double[2 * l]; 1320 s hort[] y = new short[2 * l];1165 sbyte[] y = new sbyte[2 * l]; 1321 1166 int i; 1322 1167 … … 1333 1178 1334 1179 Solver s = new Solver(); 1335 s.Solve(2 * l, new SVR_Q(prob, param), linear_term, y, 1336 alpha2, param.C, param.C, param.EPS, si, param.Shrinking); 1180 s.Solve(2 * l, new SVR_Q(prob, param), linear_term, y, alpha2, param.C, param.C, param.EPS, si, param.Shrinking); 1337 1181 1338 1182 double sum_alpha = 0; … … 1342 1186 sum_alpha += Math.Abs(alpha[i]); 1343 1187 } 1344 Debug.Write("nu = " + sum_alpha / (param.C * l) + "\n");1188 Procedures.info("nu = " + sum_alpha / (param.C * l) + "\n"); 1345 1189 } 1346 1190 … … 1352 1196 double[] alpha2 = new double[2 * l]; 1353 1197 double[] linear_term = new double[2 * l]; 1354 s hort[] y = new short[2 * l];1198 sbyte[] y = new sbyte[2 * l]; 1355 1199 int i; 1356 1200 … … 1371 1215 s.Solve(2 * l, new SVR_Q(prob, param), linear_term, y, alpha2, C, C, param.EPS, si, param.Shrinking); 1372 1216 1373 Debug.Write("epsilon = " + (-si.r) + "\n");1217 Procedures.info("epsilon = " + (-si.r) + "\n"); 1374 1218 1375 1219 for (i = 0; i < l; i++) … … 1380 1224 // decision_function 1381 1225 // 1382 privateclass decision_function1226 internal class decision_function 1383 1227 { 1384 1228 public double[] alpha; … … 1386 1230 }; 1387 1231 1388 static decision_function svm_train_one( 1389 Problem prob, Parameter param, 1390 double Cp, double Cn) 1232 static decision_function svm_train_one(Problem prob, Parameter param, double Cp, double Cn) 1391 1233 { 1392 1234 double[] alpha = new double[prob.Count]; … … 1411 1253 } 1412 1254 1413 Debug.Write("obj = " + si.obj + ", rho = " + si.rho + "\n");1255 Procedures.info("obj = " + si.obj + ", rho = " + si.rho + "\n"); 1414 1256 1415 1257 // output SVs … … 1435 1277 } 1436 1278 1437 Debug.Write("nSV = " + nSV + ", nBSV = " + nBSV + "\n");1279 Procedures.info("nSV = " + nSV + ", nBSV = " + nBSV + "\n"); 1438 1280 1439 1281 decision_function f = new decision_function(); … … 1455 1297 else prior0 += 1; 1456 1298 1457 int max_iter = 100;// Maximal number of iterations1458 double min_step = 1e-10; // Minimal step taken in line search1459 double sigma = 1e- 3; // For numerically strict PD of Hessian1299 int Max_iter = 100; // Maximal number of iterations 1300 double Min_step = 1e-10; // Minimal step taken in line search 1301 double sigma = 1e-12; // For numerically strict PD of Hessian 1460 1302 double eps = 1e-5; 1461 1303 double hiTarget = (prior1 + 1.0) / (prior1 + 2.0); … … 1480 1322 fval += (t[i] - 1) * fApB + Math.Log(1 + Math.Exp(fApB)); 1481 1323 } 1482 for (iter = 0; iter < max_iter; iter++)1324 for (iter = 0; iter < Max_iter; iter++) 1483 1325 { 1484 1326 // Update Gradient and Hessian (use H' = H + sigma I) … … 1519 1361 1520 1362 1521 stepsize = 1; 1522 while (stepsize >= min_step)1363 stepsize = 1; // Line Search 1364 while (stepsize >= Min_step) 1523 1365 { 1524 1366 newA = A + stepsize * dA; … … 1545 1387 } 1546 1388 1547 if (stepsize < min_step)1548 { 1549 Debug.Write("Line search fails in two-class probability estimates\n");1389 if (stepsize < Min_step) 1390 { 1391 Procedures.info("Line search fails in two-class probability estimates\n"); 1550 1392 break; 1551 1393 } 1552 1394 } 1553 1395 1554 if (iter >= max_iter)1555 Debug.Write("Reaching maximal iterations in two-class probability estimates\n");1396 if (iter >= Max_iter) 1397 Procedures.info("Reaching Maximal iterations in two-class probability estimates\n"); 1556 1398 probAB[0] = A; probAB[1] = B; 1557 1399 } … … 1568 1410 // Method 2 from the multiclass_prob paper by Wu, Lin, and Weng 1569 1411 private static void multiclass_probability(int k, double[,] r, double[] p) 1570 1571 int t,j;1572 int iter = 0, max_iter=Math.Max(100,k);1573 double[,] Q=new double[k,k];1574 double[] Qp= new double[k];1575 double pQp, eps=0.005/k;1576 1577 for (t=0;t<k;t++)1578 1579 p[t]=1.0/k; // Valid if k = 11580 Q[t,t]=0;1581 for (j=0;j<t;j++)1582 1583 Q[t,t]+=r[j,t]*r[j,t];1584 Q[t,j]=Q[j,t];1585 1586 for (j=t+1;j<k;j++)1587 1588 Q[t,t]+=r[j,t]*r[j,t];1589 Q[t,j]=-r[j,t]*r[t,j];1590 1591 1592 for (iter=0;iter<max_iter;iter++)1593 1594 1595 pQp=0;1596 for (t=0;t<k;t++)1597 1598 Qp[t]=0;1599 for (j=0;j<k;j++)1600 Qp[t]+=Q[t,j]*p[j];1601 pQp+=p[t]*Qp[t];1602 1603 double max_error=0;1604 for (t=0;t<k;t++)1605 1606 double error=Math.Abs(Qp[t]-pQp);1607 if (error>max_error)1608 max_error=error;1609 1610 if (max_error<eps) break;1611 1612 for (t=0;t<k;t++)1613 1614 double diff=(-Qp[t]+pQp)/Q[t,t];1615 p[t]+=diff;1616 pQp=(pQp+diff*(diff*Q[t,t]+2*Qp[t]))/(1+diff)/(1+diff);1617 for (j=0;j<k;j++)1618 1619 Qp[j]=(Qp[j]+diff*Q[t,j])/(1+diff);1620 p[j]/=(1+diff);1621 1622 1623 1624 if (iter>=max_iter)1625 Debug.Write("Exceeds max_iter in multiclass_prob\n");1626 1412 { 1413 int t, j; 1414 int iter = 0, Max_iter = Math.Max(100, k); 1415 double[,] Q = new double[k,k]; 1416 double[] Qp = new double[k]; 1417 double pQp, eps = 0.005 / k; 1418 1419 for (t = 0; t < k; t++) 1420 { 1421 p[t] = 1.0 / k; // Valid if k = 1 1422 Q[t,t] = 0; 1423 for (j = 0; j < t; j++) 1424 { 1425 Q[t,t] += r[j,t] * r[j,t]; 1426 Q[t,j] = Q[j,t]; 1427 } 1428 for (j = t + 1; j < k; j++) 1429 { 1430 Q[t,t] += r[j,t] * r[j,t]; 1431 Q[t,j] = -r[j,t] * r[t,j]; 1432 } 1433 } 1434 for (iter = 0; iter < Max_iter; iter++) 1435 { 1436 // stopping condition, recalculate QP,pQP for numerical accuracy 1437 pQp = 0; 1438 for (t = 0; t < k; t++) 1439 { 1440 Qp[t] = 0; 1441 for (j = 0; j < k; j++) 1442 Qp[t] += Q[t,j] * p[j]; 1443 pQp += p[t] * Qp[t]; 1444 } 1445 double Max_error = 0; 1446 for (t = 0; t < k; t++) 1447 { 1448 double error = Math.Abs(Qp[t] - pQp); 1449 if (error > Max_error) 1450 Max_error = error; 1451 } 1452 if (Max_error < eps) break; 1453 1454 for (t = 0; t < k; t++) 1455 { 1456 double diff = (-Qp[t] + pQp) / Q[t,t]; 1457 p[t] += diff; 1458 pQp = (pQp + diff * (diff * Q[t,t] + 2 * Qp[t])) / (1 + diff) / (1 + diff); 1459 for (j = 0; j < k; j++) 1460 { 1461 Qp[j] = (Qp[j] + diff * Q[t,j]) / (1 + diff); 1462 p[j] /= (1 + diff); 1463 } 1464 } 1465 } 1466 if (iter >= Max_iter) 1467 Procedures.info("Exceeds Max_iter in multiclass_prob\n"); 1468 } 1627 1469 1628 1470 // Cross-validation decision values for probability estimates 1629 1471 private static void svm_binary_svc_probability(Problem prob, Parameter param, double Cp, double Cn, double[] probAB) 1630 1472 { 1631 Random rand = new Random();1632 1473 int i; 1633 1474 int nr_fold = 5; … … 1636 1477 1637 1478 // random shuffle 1479 Random rand = new Random(); 1638 1480 for (i = 0; i < prob.Count; i++) perm[i] = i; 1639 1481 for (i = 0; i < prob.Count; i++) … … 1687 1529 subparam.Probability = false; 1688 1530 subparam.C = 1.0; 1689 subparam.WeightCount = 2; 1690 subparam.WeightLabels = new int[2]; 1691 subparam.Weights = new double[2]; 1692 subparam.WeightLabels[0] = +1; 1693 subparam.WeightLabels[1] = -1; 1694 subparam.Weights[0] = Cp; 1695 subparam.Weights[1] = Cn; 1531 subparam.Weights[1] = Cp; 1532 subparam.Weights[-1] = Cn; 1696 1533 Model submodel = svm_train(subprob, subparam); 1697 1534 for (j = begin; j < end; j++) … … 1718 1555 Parameter newparam = (Parameter)param.Clone(); 1719 1556 newparam.Probability = false; 1720 svm_cross_validation(prob, newparam, nr_fold, ymv , null);1557 svm_cross_validation(prob, newparam, nr_fold, ymv); 1721 1558 for (i = 0; i < prob.Count; i++) 1722 1559 { … … 1734 1571 mae += Math.Abs(ymv[i]); 1735 1572 mae /= (prob.Count - count); 1736 Debug.Write("Prob. model for test data: target value = predicted value + z,\nz: Laplace distribution e^(-|z|/sigma)/(2sigma),sigma=" + mae + "\n");1573 Procedures.info("Prob. model for test data: target value = predicted value + z,\nz: Laplace distribution e^(-|z|/sigma)/(2sigma),sigma=" + mae + "\n"); 1737 1574 return mae; 1738 1575 } … … 1743 1580 { 1744 1581 int l = prob.Count; 1745 int max_nr_class = 16;1582 int Max_nr_class = 16; 1746 1583 int nr_class = 0; 1747 int[] label = new int[ max_nr_class];1748 int[] count = new int[ max_nr_class];1584 int[] label = new int[Max_nr_class]; 1585 int[] count = new int[Max_nr_class]; 1749 1586 int[] data_label = new int[l]; 1750 1587 int i; … … 1765 1602 if (j == nr_class) 1766 1603 { 1767 if (nr_class == max_nr_class)1768 { 1769 max_nr_class *= 2;1770 int[] new_data = new int[ max_nr_class];1604 if (nr_class == Max_nr_class) 1605 { 1606 Max_nr_class *= 2; 1607 int[] new_data = new int[Max_nr_class]; 1771 1608 Array.Copy(label, 0, new_data, 0, label.Length); 1772 1609 label = new_data; 1773 new_data = new int[ max_nr_class];1610 new_data = new int[Max_nr_class]; 1774 1611 Array.Copy(count, 0, new_data, 0, count.Length); 1775 1612 count = new_data; … … 1808 1645 model.Parameter = param; 1809 1646 1810 if (param.SvmType == SvmType.ONE_CLASS || 1647 if (param.SvmType == SvmType.ONE_CLASS || 1811 1648 param.SvmType == SvmType.EPSILON_SVR || 1812 1649 param.SvmType == SvmType.NU_SVR) 1813 1650 { 1814 1651 // regression or one-class-svm 1815 model.NumberOfClasses = 2; 1652 model.NumberOfClasses = 2; 1816 1653 model.ClassLabels = null; 1817 1654 model.NumberOfSVPerClass = null; 1818 1655 model.PairwiseProbabilityA = null; model.PairwiseProbabilityB = null; 1819 1656 model.SupportVectorCoefficients = new double[1][]; 1820 1657 1821 1658 if (param.Probability && 1822 1659 (param.SvmType == SvmType.EPSILON_SVR || … … 1834 1671 int i; 1835 1672 for (i = 0; i < prob.Count; i++) 1836 if (Math.Abs(f.alpha[i]) > 0) ++nSV; 1673 if (Math.Abs(f.alpha[i]) > 0) ++nSV; 1837 1674 model.SupportVectorCount = nSV; 1838 1675 model.SupportVectors = new Node[nSV][]; … … 1873 1710 for (i = 0; i < nr_class; i++) 1874 1711 weighted_C[i] = param.C; 1875 for (i = 0; i < param.WeightCount; i++) 1876 { 1877 int j; 1878 for (j = 0; j < nr_class; j++) 1879 if (param.WeightLabels[i] == label[j]) 1880 break; 1881 if (j == nr_class) 1882 Debug.Write("warning: class label " + param.WeightLabels[i] + " specified in weight is not found\n"); 1883 else 1884 weighted_C[j] *= param.Weights[i]; 1712 foreach (int weightedLabel in param.Weights.Keys) 1713 { 1714 int index = Array.IndexOf<int>(label, weightedLabel); 1715 if (index < 0) 1716 Console.Error.WriteLine("warning: class label " + weightedLabel + " specified in weight is not found"); 1717 else weighted_C[index] *= param.Weights[weightedLabel]; 1885 1718 } 1886 1719 … … 1983 1816 } 1984 1817 1985 Debug.Write("Total nSV = " + nnz + "\n");1818 Procedures.info("Total nSV = " + nnz + "\n"); 1986 1819 1987 1820 model.SupportVectorCount = nnz; … … 2029 1862 2030 1863 // Stratified cross validation 2031 public static void svm_cross_validation(Problem prob, Parameter param, int nr_fold, double[] target , Dictionary<int,double>[] confidence)1864 public static void svm_cross_validation(Problem prob, Parameter param, int nr_fold, double[] target) 2032 1865 { 2033 1866 Random rand = new Random(); … … 2131 1964 param.SvmType == SvmType.NU_SVC)) 2132 1965 { 1966 double[] prob_estimates = new double[svm_get_nr_class(submodel)]; 2133 1967 for (j = begin; j < end; j++) 2134 {2135 double[] prob_estimates = new double[svm_get_nr_class(submodel)];2136 1968 target[perm[j]] = svm_predict_probability(submodel, prob.X[perm[j]], prob_estimates); 2137 confidence[perm[j]] = new Dictionary<int, double>();2138 for (int label = 0; label < prob_estimates.Length; label++)2139 confidence[perm[j]][submodel.ClassLabels[label]] = prob_estimates[label];2140 2141 }2142 1969 } 2143 1970 else … … 2171 1998 else 2172 1999 { 2173 Debug.Write("Model doesn't contain information for SVR probability inference\n");2000 Console.Error.WriteLine("Model doesn't contain information for SVR probability inference"); 2174 2001 return 0; 2175 2002 } … … 2185 2012 double sum = 0; 2186 2013 for (int i = 0; i < model.SupportVectorCount; i++) 2187 sum += sv_coef[i] * Kernel. k_function(x, model.SupportVectors[i], model.Parameter);2014 sum += sv_coef[i] * Kernel.KernelFunction(x, model.SupportVectors[i], model.Parameter); 2188 2015 sum -= model.Rho[0]; 2189 2016 dec_values[0] = sum; … … 2197 2024 double[] kvalue = new double[l]; 2198 2025 for (i = 0; i < l; i++) 2199 kvalue[i] = Kernel. k_function(x, model.SupportVectors[i], model.Parameter);2026 kvalue[i] = Kernel.KernelFunction(x, model.SupportVectors[i], model.Parameter); 2200 2027 2201 2028 int[] start = new int[nr_class]; … … 2262 2089 } 2263 2090 2264 int vote_ max_idx = 0;2091 int vote_Max_idx = 0; 2265 2092 for (i = 1; i < nr_class; i++) 2266 if (vote[i] > vote[vote_ max_idx])2267 vote_ max_idx = i;2268 return model.ClassLabels[vote_ max_idx];2093 if (vote[i] > vote[vote_Max_idx]) 2094 vote_Max_idx = i; 2095 return model.ClassLabels[vote_Max_idx]; 2269 2096 } 2270 2097 } 2271 2098 2272 2099 public static double svm_predict_probability(Model model, Node[] x, double[] prob_estimates) 2273 { 2274 if ((model.Parameter.SvmType == SvmType.C_SVC || model.Parameter.SvmType == SvmType.NU_SVC) && 2275 model.PairwiseProbabilityA!=null && model.PairwiseProbabilityB!=null) 2276 { 2277 int i; 2278 int nr_class = model.NumberOfClasses; 2279 double[] dec_values = new double[nr_class*(nr_class-1)/2]; 2280 svm_predict_values(model, x, dec_values); 2281 2282 double min_prob=1e-7; 2283 double[,] pairwise_prob=new double[nr_class,nr_class]; 2284 2285 int k=0; 2286 for(i=0;i<nr_class;i++) 2287 for(int j=i+1;j<nr_class;j++) 2288 { 2289 pairwise_prob[i,j]=Math.Min(Math.Max(sigmoid_predict(dec_values[k],model.PairwiseProbabilityA[k],model.PairwiseProbabilityB[k]),min_prob),1-min_prob); 2290 pairwise_prob[j,i]=1-pairwise_prob[i,j]; 2291 k++; 2292 } 2293 multiclass_probability(nr_class,pairwise_prob,prob_estimates); 2294 2295 int prob_max_idx = 0; 2296 for(i=1;i<nr_class;i++) 2297 if(prob_estimates[i] > prob_estimates[prob_max_idx]) 2298 prob_max_idx = i; 2299 return model.ClassLabels[prob_max_idx]; 2300 } 2301 else 2302 return svm_predict(model, x); 2303 } 2304 2305 private static double atof(string s) 2306 { 2307 return double.Parse(s); 2308 } 2309 2310 private static int atoi(string s) 2311 { 2312 return int.Parse(s); 2100 { 2101 if ((model.Parameter.SvmType == SvmType.C_SVC || model.Parameter.SvmType == SvmType.NU_SVC) && 2102 model.PairwiseProbabilityA != null && model.PairwiseProbabilityB != null) 2103 { 2104 int i; 2105 int nr_class = model.NumberOfClasses; 2106 double[] dec_values = new double[nr_class * (nr_class - 1) / 2]; 2107 svm_predict_values(model, x, dec_values); 2108 2109 double Min_prob = 1e-7; 2110 double[,] pairwise_prob = new double[nr_class, nr_class]; 2111 2112 int k = 0; 2113 for (i = 0; i < nr_class; i++) 2114 { 2115 for (int j = i + 1; j < nr_class; j++) 2116 { 2117 pairwise_prob[i, j] = Math.Min(Math.Max(sigmoid_predict(dec_values[k], model.PairwiseProbabilityA[k], model.PairwiseProbabilityB[k]), Min_prob), 1 - Min_prob); 2118 pairwise_prob[j, i] = 1 - pairwise_prob[i, j]; 2119 k++; 2120 } 2121 } 2122 multiclass_probability(nr_class, pairwise_prob, prob_estimates); 2123 2124 int prob_Max_idx = 0; 2125 for (i = 1; i < nr_class; i++) 2126 if (prob_estimates[i] > prob_estimates[prob_Max_idx]) 2127 prob_Max_idx = i; 2128 return model.ClassLabels[prob_Max_idx]; 2129 } 2130 else 2131 return svm_predict(model, x); 2313 2132 } 2314 2133 … … 2318 2137 2319 2138 SvmType svm_type = param.SvmType; 2320 if (svm_type != SvmType.C_SVC &&2321 svm_type != SvmType.NU_SVC &&2322 svm_type != SvmType.ONE_CLASS &&2323 svm_type != SvmType.EPSILON_SVR &&2324 svm_type != SvmType.NU_SVR)2325 return "unknown svm type";2326 2139 2327 2140 // kernel_type, degree 2328 2141 2329 2142 KernelType kernel_type = param.KernelType; 2330 if (kernel_type != KernelType.LINEAR &&2331 kernel_type != KernelType.POLY &&2332 kernel_type != KernelType.RBF &&2333 kernel_type != KernelType.SIGMOID &&2334 kernel_type != KernelType.PRECOMPUTED)2335 return "unknown kernel type";2336 2143 2337 2144 if (param.Degree < 0) … … 2345 2152 if (param.EPS <= 0) 2346 2153 return "eps <= 0"; 2154 2155 if (param.Gamma == 0) 2156 param.Gamma = 1.0 / prob.MaxIndex; 2347 2157 2348 2158 if (svm_type == SvmType.C_SVC || … … 2362 2172 return "p < 0"; 2363 2173 2364 if (param.Probability && svm_type == SvmType.ONE_CLASS) 2174 if (param.Probability && 2175 svm_type == SvmType.ONE_CLASS) 2365 2176 return "one-class SVM probability output not supported yet"; 2366 2177 … … 2370 2181 { 2371 2182 int l = prob.Count; 2372 int max_nr_class = 16;2183 int Max_nr_class = 16; 2373 2184 int nr_class = 0; 2374 int[] label = new int[ max_nr_class];2375 int[] count = new int[ max_nr_class];2185 int[] label = new int[Max_nr_class]; 2186 int[] count = new int[Max_nr_class]; 2376 2187 2377 2188 int i; … … 2389 2200 if (j == nr_class) 2390 2201 { 2391 if (nr_class == max_nr_class)2392 { 2393 max_nr_class *= 2;2394 int[] new_data = new int[ max_nr_class];2202 if (nr_class == Max_nr_class) 2203 { 2204 Max_nr_class *= 2; 2205 int[] new_data = new int[Max_nr_class]; 2395 2206 Array.Copy(label, 0, new_data, 0, label.Length); 2396 2207 label = new_data; 2397 2208 2398 new_data = new int[ max_nr_class];2209 new_data = new int[Max_nr_class]; 2399 2210 Array.Copy(count, 0, new_data, 0, count.Length); 2400 2211 count = new_data; … … 2432 2243 } 2433 2244 } 2434 2435 2245 } -
trunk/sources/LibSVM/Training.cs
r1819 r2415 23 23 namespace SVM 24 24 { 25 /// < remarks>25 /// <summary> 26 26 /// Class containing the routines to train SVM models. 27 /// </ remarks>27 /// </summary> 28 28 public static class Training 29 29 { 30 /// <summary> 31 /// Whether the system will output information to the console during the training process. 32 /// </summary> 33 public static bool IsVerbose 34 { 35 get 36 { 37 return Procedures.IsVerbose; 38 } 39 set 40 { 41 Procedures.IsVerbose = value; 42 } 43 } 44 30 45 private static double doCrossValidation(Problem problem, Parameter parameters, int nr_fold) 31 46 { 32 47 int i; 33 48 double[] target = new double[problem.Count]; 34 Dictionary<int, double>[] confidence = new Dictionary<int, double>[problem.Count]; 35 Procedures.svm_cross_validation(problem, parameters, nr_fold, target, confidence); 36 if (parameters.Probability) 37 { 38 List<RankPair> ranks = new List<RankPair>(); 39 for (i = 0; i < target.Length; i++) 40 ranks.Add(new RankPair(confidence[i][1], problem.Y[i])); 41 PerformanceEvaluator eval = new PerformanceEvaluator(ranks); 42 return eval.AuC*eval.AP; 49 Procedures.svm_cross_validation(problem, parameters, nr_fold, target); 50 int total_correct = 0; 51 double total_error = 0; 52 double sumv = 0, sumy = 0, sumvv = 0, sumyy = 0, sumvy = 0; 53 if (parameters.SvmType == SvmType.EPSILON_SVR || parameters.SvmType == SvmType.NU_SVR) 54 { 55 for (i = 0; i < problem.Count; i++) 56 { 57 double y = problem.Y[i]; 58 double v = target[i]; 59 total_error += (v - y) * (v - y); 60 sumv += v; 61 sumy += y; 62 sumvv += v * v; 63 sumyy += y * y; 64 sumvy += v * y; 65 } 66 return(problem.Count * sumvy - sumv * sumy) / (Math.Sqrt(problem.Count * sumvv - sumv * sumv) * Math.Sqrt(problem.Count * sumyy - sumy * sumy)); 43 67 } 44 68 else 45 { 46 int total_correct = 0; 47 double total_error = 0; 48 double sumv = 0, sumy = 0, sumvv = 0, sumyy = 0, sumvy = 0; 49 if (parameters.SvmType == SvmType.EPSILON_SVR || parameters.SvmType == SvmType.NU_SVR) 50 { 51 for (i = 0; i < problem.Count; i++) 52 { 53 double y = problem.Y[i]; 54 double v = target[i]; 55 total_error += (v - y) * (v - y); 56 sumv += v; 57 sumy += y; 58 sumvv += v * v; 59 sumyy += y * y; 60 sumvy += v * y; 61 } 62 } 63 else 64 for (i = 0; i < problem.Count; i++) 65 if (target[i] == problem.Y[i]) 66 ++total_correct; 67 return (double)total_correct / problem.Count; 68 } 69 69 for (i = 0; i < problem.Count; i++) 70 if (target[i] == problem.Y[i]) 71 ++total_correct; 72 return (double)total_correct / problem.Count; 70 73 } 71 74 /// <summary> … … 96 99 public static double PerformCrossValidation(Problem problem, Parameter parameters, int nrfold) 97 100 { 98 Procedures.svm_check_parameter(problem, parameters); 99 return doCrossValidation(problem, parameters, nrfold); 101 string error = Procedures.svm_check_parameter(problem, parameters); 102 if (error == null) 103 return doCrossValidation(problem, parameters, nrfold); 104 else throw new Exception(error); 100 105 } 101 106 … … 108 113 public static Model Train(Problem problem, Parameter parameters) 109 114 { 110 Procedures.svm_check_parameter(problem, parameters); 111 112 return Procedures.svm_train(problem, parameters); 115 string error = Procedures.svm_check_parameter(problem, parameters); 116 117 if (error == null) 118 return Procedures.svm_train(problem, parameters); 119 else throw new Exception(error); 113 120 } 114 121 … … 190 197 191 198 case 'w': 192 ++parameters.WeightCount; 193 { 194 int[] old = parameters.WeightLabels; 195 parameters.WeightLabels = new int[parameters.WeightCount]; 196 Array.Copy(old, 0, parameters.WeightLabels, 0, parameters.WeightCount - 1); 197 } 198 { 199 double[] old = parameters.Weights; 200 parameters.Weights = new double[parameters.WeightCount]; 201 Array.Copy(old, 0, parameters.Weights, 0, parameters.WeightCount - 1); 202 } 203 204 parameters.WeightLabels[parameters.WeightCount - 1] = int.Parse(args[i - 1].Substring(2)); 205 parameters.Weights[parameters.WeightCount - 1] = double.Parse(args[i]); 199 parameters.Weights[int.Parse(args[i - 1].Substring(2))] = double.Parse(args[1]); 206 200 break; 207 201
Note: See TracChangeset
for help on using the changeset viewer.