4 files changed, 293 insertions, 73 deletions
diff --git a/lib/classifier/lsi/content_node.rb b/lib/classifier/lsi/content_node.rb
index f313331..37a9448 100644
--- a/lib/classifier/lsi/content_node.rb
+++ b/lib/classifier/lsi/content_node.rb
@@ -1,72 +1,96 @@
+# rbs_inline: enabled
+
 # Author::    David Fayram  (mailto:dfayram@lensmen.net)
 # Copyright:: Copyright (c) 2005 David Fayram II
 # License::   LGPL
 
 module Classifier
-
-# This is an internal data structure class for the LSI node. Save for 
-# raw_vector_with, it should be fairly straightforward to understand.
-# You should never have to use it directly.
+  # This is an internal data structure class for the LSI node. Save for
+  # raw_vector_with, it should be fairly straightforward to understand.
+  # You should never have to use it directly.
   class ContentNode
-    attr_accessor :raw_vector, :raw_norm, 
-                  :lsi_vector, :lsi_norm,
-                  :categories   
-    
+    # @rbs @word_hash: Hash[Symbol, Integer]
+
+    # @rbs @raw_vector: untyped
+    # @rbs @raw_norm: untyped
+    # @rbs @lsi_vector: untyped
+    # @rbs @lsi_norm: untyped
+    attr_accessor :raw_vector, :raw_norm, :lsi_vector, :lsi_norm
+
+    # @rbs @categories: Array[String | Symbol]
+    attr_accessor :categories
+
     attr_reader :word_hash
+
     # If text_proc is not specified, the source will be duck-typed
     # via source.to_s
-    def initialize( word_hash, *categories )
+    #
+    # @rbs (Hash[Symbol, Integer], *String | Symbol) -> void
+    def initialize(word_frequencies, *categories)
       @categories = categories || []
-      @word_hash = word_hash
+      @word_hash = word_frequencies
     end
-   
+
     # Use this to fetch the appropriate search vector.
+    #
+    # @rbs () -> untyped
     def search_vector
       @lsi_vector || @raw_vector
     end
-    
+
     # Use this to fetch the appropriate search vector in normalized form.
+    #
+    # @rbs () -> untyped
     def search_norm
       @lsi_norm || @raw_norm
     end
-   
+
     # Creates the raw vector out of word_hash using word_list as the
     # key for mapping the vector space.
-    def raw_vector_with( word_list )
-      if $GSL
-         vec = GSL::Vector.alloc(word_list.size)
-      else
-         vec = Array.new(word_list.size, 0)
-      end
+    #
+    # @rbs (WordList) -> untyped
+    def raw_vector_with(word_list)
+      vec = if Classifier::LSI.native_available?
+              Classifier::LSI.vector_class.alloc(word_list.size)
+            else
+              Array.new(word_list.size, 0)
+            end
 
       @word_hash.each_key do |word|
         vec[word_list[word]] = @word_hash[word] if word_list[word]
       end
-     
+
       # Perform the scaling transform
-      total_words = vec.sum
-      
+      total_words = Classifier::LSI.native_available? ? vec.sum : vec.sum_with_identity
+      vec_array = Classifier::LSI.native_available? ? vec.to_a : vec
+      total_unique_words = vec_array.count { |word| word != 0 }
+
       # Perform first-order association transform if this vector has more
-      # than one word in it. 
-      if total_words > 1.0 
+      # than one word in it.
+      if total_words > 1.0 && total_unique_words > 1
         weighted_total = 0.0
+
         vec.each do |term|
-          if ( term > 0 )
-            weighted_total += (( term / total_words ) * Math.log( term / total_words ))
-          end
-        end 
-        vec = vec.collect { |val| Math.log( val + 1 ) / -weighted_total }
+          next unless term.positive?
+          next if total_words.zero?
+
+          term_over_total = term / total_words
+          val = term_over_total * Math.log(term_over_total)
+          weighted_total += val unless val.nan?
+        end
+
+        sign = weighted_total.negative? ? 1.0 : -1.0
+        divisor = sign * [weighted_total.abs, Vector::EPSILON].max
+        vec = vec.collect { |val| Math.log(val + 1) / divisor }
       end
-      
-      if $GSL
-         @raw_norm   = vec.normalize
-         @raw_vector = vec
+
+      if Classifier::LSI.native_available?
+        @raw_norm   = vec.normalize
+        @raw_vector = vec
       else
-         @raw_norm   = Vector[*vec].normalize
-         @raw_vector = Vector[*vec]
+        @raw_norm   = Vector[*vec].normalize
+        @raw_vector = Vector[*vec]
       end
-    end   
-  
-  end  
-  
+    end
+  end
 end
diff --git a/lib/classifier/lsi/incremental_svd.rb b/lib/classifier/lsi/incremental_svd.rb
new file mode 100644
index 0000000..beabfa3
--- /dev/null
+++ b/lib/classifier/lsi/incremental_svd.rb
@@ -0,0 +1,166 @@
+# rbs_inline: enabled
+
+# rubocop:disable Naming/MethodParameterName, Metrics/ParameterLists
+
+require 'matrix'
+
+module Classifier
+  class LSI
+    # Brand's Incremental SVD Algorithm for LSI
+    #
+    # Implements the algorithm from Brand (2006) "Fast low-rank modifications
+    # of the thin singular value decomposition" for adding documents to LSI
+    # without full SVD recomputation.
+    #
+    # Given existing thin SVD: A ≈ U * S * V^T (with k components)
+    # When adding a new column c:
+    #
+    # 1. Project: m = U^T * c (project onto existing column space)
+    # 2. Residual: p = c - U * m (component orthogonal to U)
+    # 3. Orthonormalize: If ||p|| > ε: p̂ = p / ||p||
+    # 4. Form K matrix:
+    #    - If ||p|| > ε: K = [diag(s), m; 0, ||p||] (rank grows by 1)
+    #    - If ||p|| ≈ 0: K = diag(s) + m * e_last^T (no new direction)
+    # 5. Small SVD: Compute SVD of K (only (k+1) × (k+1) matrix!)
+    # 6. Update:
+    #    - U_new = [U, p̂] * U'
+    #    - S_new = S'
+    #
+    module IncrementalSVD
+      EPSILON = 1e-10
+
+      class << self
+        # Updates SVD with a new document vector using Brand's algorithm.
+        #
+        # @param u [Matrix] current left singular vectors (m × k)
+        # @param s [Array<Float>] current singular values (k values)
+        # @param c [Vector] new document vector (m × 1)
+        # @param max_rank [Integer] maximum rank to maintain
+        # @param epsilon [Float] threshold for zero detection
+        # @return [Array<Matrix, Array<Float>>] updated [u, s]
+        #
+        # @rbs (Matrix, Array[Float], Vector, max_rank: Integer, ?epsilon: Float) -> [Matrix, Array[Float]]
+        def update(u, s, c, max_rank:, epsilon: EPSILON)
+          m_vec = project(u, c)
+          u_times_m = u * m_vec
+          p_vec = c - (u_times_m.is_a?(Vector) ? u_times_m : Vector[*u_times_m.to_a.flatten])
+          p_norm = magnitude(p_vec)
+
+          if p_norm > epsilon
+            update_with_new_direction(u, s, m_vec, p_vec, p_norm, max_rank, epsilon)
+          else
+            update_in_span(u, s, m_vec, max_rank, epsilon)
+          end
+        end
+
+        # Projects a document vector onto the semantic space defined by U.
+        # Returns the LSI representation: lsi_vec = U^T * raw_vec
+        #
+        # @param u [Matrix] left singular vectors (m × k)
+        # @param raw_vec [Vector] document vector in term space (m × 1)
+        # @return [Vector] document in semantic space (k × 1)
+        #
+        # @rbs (Matrix, Vector) -> Vector
+        def project(u, raw_vec)
+          u.transpose * raw_vec
+        end
+
+        private
+
+        # Update when new document has a component orthogonal to existing U.
+        # @rbs (Matrix, Array[Float], Vector, Vector, Float, Integer, Float) -> [Matrix, Array[Float]]
+        def update_with_new_direction(u, s, m_vec, p_vec, p_norm, max_rank, epsilon)
+          p_hat = p_vec * (1.0 / p_norm)
+          k_matrix = build_k_matrix_with_growth(s, m_vec, p_norm)
+          u_prime, s_prime = small_svd(k_matrix, epsilon)
+          u_extended = extend_matrix_with_column(u, p_hat)
+          u_new = u_extended * u_prime
+
+          u_new, s_prime = truncate(u_new, s_prime, max_rank) if s_prime.size > max_rank
+
+          [u_new, s_prime]
+        end
+
+        # Update when new document is entirely within span of existing U.
+        # @rbs (Matrix, Array[Float], Vector, Integer, Float) -> [Matrix, Array[Float]]
+        def update_in_span(u, s, m_vec, max_rank, epsilon)
+          k_matrix = build_k_matrix_in_span(s, m_vec)
+          u_prime, s_prime = small_svd(k_matrix, epsilon)
+          u_new = u * u_prime
+
+          u_new, s_prime = truncate(u_new, s_prime, max_rank) if s_prime.size > max_rank
+
+          [u_new, s_prime]
+        end
+
+        # Builds the K matrix when rank grows by 1.
+        # @rbs (Array[Float], untyped, Float) -> untyped
+        def build_k_matrix_with_growth(s, m_vec, p_norm)
+          k = s.size
+          rows = k.times.map do |i|
+            row = Array.new(k + 1, 0.0) #: Array[Float]
+            row[i] = s[i].to_f
+            row[k] = m_vec[i].to_f
+            row
+          end
+          rows << Array.new(k + 1, 0.0).tap { |r| r[k] = p_norm }
+          Matrix.rows(rows)
+        end
+
+        # Builds the K matrix when vector is in span (no rank growth).
+        # @rbs (Array[Float], Vector) -> Matrix
+        def build_k_matrix_in_span(s, _m_vec)
+          k = s.size
+          rows = k.times.map do |i|
+            row = Array.new(k, 0.0)
+            row[i] = s[i]
+            row
+          end
+          Matrix.rows(rows)
+        end
+
+        # Computes SVD of small matrix and extracts singular values.
+        # @rbs (Matrix, Float) -> [Matrix, Array[Float]]
+        def small_svd(matrix, epsilon)
+          u, _v, s_array = matrix.SV_decomp
+
+          s_sorted = s_array.select { |sv| sv.abs > epsilon }.sort.reverse
+          indices = s_array.each_with_index
+                           .select { |sv, _| sv.abs > epsilon }
+                           .sort_by { |sv, _| -sv }
+                           .map { |_, i| i }
+
+          u_cols = indices.map { |i| u.column(i).to_a }
+          u_reordered = u_cols.empty? ? Matrix.empty(matrix.row_size, 0) : Matrix.columns(u_cols)
+
+          [u_reordered, s_sorted]
+        end
+
+        # Extends matrix with a new column
+        # @rbs (Matrix, Vector) -> Matrix
+        def extend_matrix_with_column(matrix, col_vec)
+          rows = matrix.row_size.times.map do |i|
+            matrix.row(i).to_a + [col_vec[i]]
+          end
+          Matrix.rows(rows)
+        end
+
+        # Truncates to max_rank
+        # @rbs (untyped, Array[Float], Integer) -> [untyped, Array[Float]]
+        def truncate(u, s, max_rank)
+          s_truncated = s[0...max_rank] || [] #: Array[Float]
+          cols = (0...max_rank).map { |i| u.column(i).to_a }
+          u_truncated = Matrix.columns(cols)
+          [u_truncated, s_truncated]
+        end
+
+        # Computes magnitude of a vector
+        # @rbs (untyped) -> Float
+        def magnitude(vec)
+          Math.sqrt(vec.to_a.sum { |x| x.to_f * x.to_f })
+        end
+      end
+    end
+  end
+end
+# rubocop:enable Naming/MethodParameterName, Metrics/ParameterLists
diff --git a/lib/classifier/lsi/summary.rb b/lib/classifier/lsi/summary.rb
index 6adf7ac..3328057 100644
--- a/lib/classifier/lsi/summary.rb
+++ b/lib/classifier/lsi/summary.rb
@@ -3,29 +3,49 @@
 # License::   LGPL
 
 class String
-   def summary( count=10, separator=" [...] " )
-      perform_lsi split_sentences, count, separator
-   end
-
-   def paragraph_summary( count=1, separator=" [...] " )
-      perform_lsi split_paragraphs, count, separator
-   end
-
-   def split_sentences
-      split /(\.|\!|\?)/ # TODO: make this less primitive
-   end
-   
-   def split_paragraphs
-      split /(\n\n|\r\r|\r\n\r\n)/ # TODO: make this less primitive
-   end
-   
-   private
-   
-   def perform_lsi(chunks, count, separator)
-      lsi = Classifier::LSI.new :auto_rebuild => false
-      chunks.each { |chunk| lsi << chunk unless chunk.strip.empty? || chunk.strip.split.size == 1 }
-      lsi.build_index
-      summaries = lsi.highest_relative_content count
-      return summaries.reject { |chunk| !summaries.include? chunk }.map { |x| x.strip }.join(separator)
-   end
-end
-\ No newline at end of file
+  ABBREVIATIONS = %w[Mr Mrs Ms Dr Prof Jr Sr Inc Ltd Corp Co vs etc al eg ie].freeze
+
+  def summary(count = 10, separator = ' [...] ')
+    perform_lsi split_sentences, count, separator
+  end
+
+  def paragraph_summary(count = 1, separator = ' [...] ')
+    perform_lsi split_paragraphs, count, separator
+  end
+
+  def split_sentences
+    return pragmatic_segment if defined?(PragmaticSegmenter)
+
+    split_sentences_regex
+  end
+
+  def split_paragraphs
+    split(/\r?\n\r?\n+/)
+  end
+
+  private
+
+  def pragmatic_segment
+    PragmaticSegmenter::Segmenter.new(text: self).segment
+  end
+
+  def split_sentences_regex
+    abbrev_pattern = ABBREVIATIONS.map { |a| "#{a}\\." }.join('|')
+    text = gsub(/\b(#{abbrev_pattern})/i) { |m| m.gsub('.', '<<<DOT>>>') }
+    text = text.gsub(/(\d)\.(\d)/, '\1<<<DOT>>>\2')
+    sentences = text.split(/(?<=[.!?])(?:\s+|(?=[A-Z]))/)
+    sentences.map { |s| s.gsub('<<<DOT>>>', '.') }
+  end
+
+  def perform_lsi(chunks, count, separator)
+    lsi = Classifier::LSI.new auto_rebuild: false
+    chunks.each do |chunk|
+      stripped = chunk.strip
+      next if stripped.empty? || stripped.split.size == 1
+
+      lsi << chunk
+    end
+    lsi.build_index
+    lsi.highest_relative_content(count).map(&:strip).join(separator)
+  end
+end
diff --git a/lib/classifier/lsi/word_list.rb b/lib/classifier/lsi/word_list.rb
index dba3bde..a8bca9e 100644
--- a/lib/classifier/lsi/word_list.rb
+++ b/lib/classifier/lsi/word_list.rb
@@ -1,36 +1,46 @@
+# rbs_inline: enabled
+
 # Author::    David Fayram  (mailto:dfayram@lensmen.net)
 # Copyright:: Copyright (c) 2005 David Fayram II
 # License::   LGPL
 
-module Classifier  
+module Classifier
   # This class keeps a word => index mapping. It is used to map stemmed words
   # to dimensions of a vector.
-    
   class WordList
+    # @rbs @location_table: Hash[Symbol, Integer]
+
+    # @rbs () -> void
     def initialize
-      @location_table = Hash.new
+      @location_table = {}
     end
-    
+
     # Adds a word (if it is new) and assigns it a unique dimension.
+    #
+    # @rbs (Symbol) -> Integer?
     def add_word(word)
       term = word
       @location_table[term] = @location_table.size unless @location_table[term]
     end
-    
+
     # Returns the dimension of the word or nil if the word is not in the space.
+    #
+    # @rbs (Symbol) -> Integer?
     def [](lookup)
       term = lookup
       @location_table[term]
     end
-    
+
+    # @rbs (Integer) -> Symbol?
     def word_for_index(ind)
       @location_table.invert[ind]
     end
-    
+
     # Returns the number of words mapped.
+    #
+    # @rbs () -> Integer
     def size
       @location_table.size
     end
-    
   end
 end