Merge pull request SciRuby#31 from v0dro/wo-gsl

changed stuff to ensure proper working with or without GSL.

lib/statsample.rb

@@ -141,7 +141,7 @@ def self.create_has_library(library)
       class_variable_get(cv)
     end
   end
-
+  
   create_has_library :gsl
 
   SPLIT_TOKEN = ','

lib/statsample/matrix.rb

@@ -10,6 +10,7 @@ class ::Matrix
   def to_matrix
     self
   end
+
   def to_dataset
     f = (self.respond_to? :fields_y) ? fields_y : column_size.times.map {|i| _("VAR_%d") % (i+1) }
     ds=Statsample::Dataset.new(f)
@@ -24,6 +25,7 @@ def to_dataset
     ds.name=self.name if self.respond_to? :name
     ds
   end
+
   if defined? :eigenpairs
     alias_method :eigenpairs_ruby, :eigenpairs
   end
@@ -38,17 +40,15 @@ def eigenpairs
   def eigenvalues
     eigenpairs.collect {|v| v[0]}
   end
+
   def eigenvectors
     eigenpairs.collect {|v| v[1]}
   end
+
   def eigenvectors_matrix
     Matrix.columns(eigenvectors)
   end
 
-
-
-
-
   def to_gsl
     out=[]
     self.row_size.times{|i|
@@ -64,9 +64,11 @@ class Col
       def to_matrix
       ::Matrix.columns([self.size.times.map {|i| self[i]}])
       end
+
       def to_ary
         to_a
       end
+
       def to_gsl
         self
       end
@@ -95,18 +97,23 @@ def to_dataset
     def row_size
       size1
     end
+
     def column_size
       size2
     end
+
     def determinant
       det
     end
+
     def inverse
       GSL::Linalg::LU.invert(self)
     end
+
     def eigenvalues
       eigenpairs.collect {|v| v[0]}
     end
+
     def eigenvectors
       eigenpairs.collect {|v| v[1]}
     end
@@ -123,6 +130,7 @@ def eigenvectors_matrix
       GSL::Eigen::symmv_sort(eigval, eigvec, GSL::Eigen::SORT_VAL_DESC)
       eigvec
     end
+
     def eigenpairs
       eigval, eigvec= GSL::Eigen.symmv(self)
       GSL::Eigen::symmv_sort(eigval, eigvec, GSL::Eigen::SORT_VAL_DESC)
@@ -137,12 +145,14 @@ def eigenpairs
     def square?
       size1==size2
     end
+
     def to_matrix
       rows=self.size1
       cols=self.size2
       out=(0...rows).collect{|i| (0...cols).collect {|j| self[i,j]} }
       ::Matrix.rows(out)
     end
+
     def total_sum
       sum=0
       size1.times {|i|

lib/statsample/test/bartlettsphericity.rb

@@ -31,7 +31,7 @@ def initialize(matrix,ncases)
       # 
       def compute
         @value=-((@ncases-1)-(2*@nvars+5).quo(6))*Math::log(@matrix.determinant)
-        @df=(@nvars*(@nvars-1)).quo(2)
+        @df=(@nvars*(@nvars-1)) / 2
       end
       def probability
         1-Distribution::ChiSquare.cdf(@value,@df)

statsample.gemspec

@@ -85,4 +85,5 @@ Gem::Specification.new do |s|
   s.add_development_dependency 'minitest', '~> 5.5'
   s.add_development_dependency 'gettext', '~> 3.1'
   s.add_development_dependency 'mocha', '~> 1.1'
+  # s.add_development_dependency 'gsl-nmatrix', '~> 1.17'
 end

test/test_factor.rb

@@ -34,51 +34,53 @@ def test_covariance_matrix
   end
 
   def test_principalcomponents_ruby_gsl
-    ran = Distribution::Normal.rng
-
-    #    @r=::Rserve::Connection.new
-
-    samples = 20
-    [3, 5, 7].each {|k|
-      v = {}
-      v['x0'] = samples.times.map { ran.call }.to_scale.centered
-      (1...k).each {|i|
-        v["x#{i}"] = samples.times.map { |ii| ran.call * 0.5 + v["x#{i - 1}"][ii] * 0.5 }.to_scale.centered
-      }
-
-      ds = v.to_dataset
-      cm = ds.covariance_matrix
-      #      @r.assign('ds',ds)
-      #      @r.eval('cm<-cor(ds);sm<-eigen(cm, sym=TRUE);v<-sm$vectors')
-      #      puts "eigenvalues"
-      #      puts @r.eval('v').to_ruby.to_s
-      pca_ruby = Statsample::Factor::PCA.new(cm, m: k, use_gsl: false)
-      pca_gsl = Statsample::Factor::PCA.new(cm, m: k, use_gsl: true)
-      pc_ruby = pca_ruby.principal_components(ds)
-      pc_gsl  = pca_gsl.principal_components(ds)
-      # Test component matrix correlation!
-      cm_ruby = pca_ruby.component_matrix
-      # puts cm_ruby.summary
-      k.times {|i|
-        pc_id = "PC_#{i + 1}"
-        assert_in_delta(pca_ruby.eigenvalues[i], pca_gsl.eigenvalues[i], 1e-10)
-        # Revert gsl component values
-        pc_gsl_data = (pc_gsl[pc_id][0] - pc_ruby[pc_id][0]).abs > 1e-6 ? pc_gsl[pc_id].recode(&:-@) : pc_gsl[pc_id]
-        assert_similar_vector(pc_gsl_data, pc_ruby[pc_id], 1e-6, "PC for #{k} variables")
-        if false
-          k.times {|j| # variable
-            ds_id = "x#{j}"
-            r = Statsample::Bivariate.correlation(ds[ds_id], pc_ruby[pc_id])
-            puts "#{pc_id}-#{ds_id}:#{r}"
-          }
-        end
+    if Statsample.has_gsl?
+      ran = Distribution::Normal.rng
+
+      #    @r=::Rserve::Connection.new
+
+      samples = 20
+      [3, 5, 7].each {|k|
+        v = {}
+        v['x0'] = samples.times.map { ran.call }.to_scale.centered
+        (1...k).each {|i|
+          v["x#{i}"] = samples.times.map { |ii| ran.call * 0.5 + v["x#{i - 1}"][ii] * 0.5 }.to_scale.centered
+        }
+
+        ds = v.to_dataset
+        cm = ds.covariance_matrix
+        #      @r.assign('ds',ds)
+        #      @r.eval('cm<-cor(ds);sm<-eigen(cm, sym=TRUE);v<-sm$vectors')
+        #      puts "eigenvalues"
+        #      puts @r.eval('v').to_ruby.to_s
+        pca_ruby = Statsample::Factor::PCA.new(cm, m: k, use_gsl: false)
+        pca_gsl = Statsample::Factor::PCA.new(cm, m: k, use_gsl: true)
+        pc_ruby = pca_ruby.principal_components(ds)
+        pc_gsl  = pca_gsl.principal_components(ds)
+        # Test component matrix correlation!
+        cm_ruby = pca_ruby.component_matrix
+        # puts cm_ruby.summary
+        k.times {|i|
+          pc_id = "PC_#{i + 1}"
+          assert_in_delta(pca_ruby.eigenvalues[i], pca_gsl.eigenvalues[i], 1e-10)
+          # Revert gsl component values
+          pc_gsl_data = (pc_gsl[pc_id][0] - pc_ruby[pc_id][0]).abs > 1e-6 ? pc_gsl[pc_id].recode(&:-@) : pc_gsl[pc_id]
+          assert_similar_vector(pc_gsl_data, pc_ruby[pc_id], 1e-6, "PC for #{k} variables")
+          if false
+            k.times {|j| # variable
+              ds_id = "x#{j}"
+              r = Statsample::Bivariate.correlation(ds[ds_id], pc_ruby[pc_id])
+              puts "#{pc_id}-#{ds_id}:#{r}"
+            }
+          end
+        }
       }
-    }
+    end
     # @r.close
   end
 
   def test_principalcomponents
-    principalcomponents(true)
+    principalcomponents(true) if Statsample.has_gsl?
     principalcomponents(false)
   end