diff --git a/README.md b/README.md index 5d244033..62d2fff2 100644 --- a/README.md +++ b/README.md @@ -187,61 +187,70 @@ puts a.matmul(a) # [15, 22]] ``` -### DataFrames +### Machine Learning -For more structured data, consider using a `DataFrame` +`Num::Grad` provides a pure-crystal approach to find derivatives of +mathematical functions. Use a `Num::Grad::Variable` with a `Num::Grad::Context` +to easily compute these derivatives. ```crystal -df = DataFrame.from_items( - foo: [1, 2, 3, 4, 5].to_tensor, - bar: [2.73, 3.1, 4.8, 5.1, 3.2], -) - -puts df - -# foo bar -# 0 1 2.73 -# 1 2 3.1 -# 2 3 4.8 -# 3 4 5.1 -# 4 5 3.2 +ctx = Num::Grad::Context(Tensor(Float64)).new + +x = ctx.variable([3.0]) +y = ctx.variable([2.0]) + +# f(x) = x ** y +f = x ** y +puts f # => [9] + +f.backprop + +# df/dx = y * x = 6.0 +puts x.grad # => [6.0] ``` -A `DataFrame` maintains types while still providing convenient -mapping and reduction operations +`Num::NN` contains an extension to `Num::Grad` that provides an easy-to-use +interface to assist in creating neural networks. Designing and creating +a network is simple using Crystal's block syntax. ```crystal -puts df.c[:foo] +ctx = Num::Grad::Context(Tensor(Float64)).new -# 0 1 -# 1 2 -# 2 3 -# 3 4 -# 4 5 -# Name: foo -# dtype: Int32 +x_train = [[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]].to_tensor +y_train = [[0.0], [1.0], [1.0], [0.0]].to_tensor -puts typeof(df.c[:foo]) +x = ctx.variable(x_train) -# Series(Int32, Int32) +net = Num::NN::Network.new(ctx) do -puts df.sum + # A basic network with a single hidden layer using + # a ReLU activation function + linear(2, 3) + relu + linear(3, 1) -# foo 15 -# bar 18.93 -``` + # SGD Optimizer + sgd 0.7 -With operations that broadcast across the `DataFrame` + # Sigmoid Cross Entropy to calculate loss + sigmoid_cross_entropy_loss +end -```crystal -puts df.greater(df.mean) - -# foo bar -# 0 false false -# 1 false false -# 2 false true -# 3 true true -# 4 true false +500.times do |epoch| + y_pred = net.forward(x) + loss = net.loss(y_pred, y_train) + puts "Epoch: #{epoch} - Loss #{loss}" + loss.backprop + net.optimizer.update +end + +# Clip results to make a prediction +puts net.forward(x).value.map { |el| el > 0 ? 1 : 0} + +# [[0], +# [1], +# [1], +# [0]] ``` Review the documentation for full implementation details, and if something is missing, diff --git a/examples/basic_xor_classifier/README.md b/examples/basic_xor_classifier/README.md new file mode 100644 index 00000000..6f4a91ba --- /dev/null +++ b/examples/basic_xor_classifier/README.md @@ -0,0 +1,74 @@ +## Basic XOR Classifier + +The following implements a simple XOR classifier to show how to use +`num.cr`'s `Network` class. Plotting is done via `ishi`. + +```crystal +ctx = Num::Grad::Context(Tensor(Float64)).new + +bsz = 32 + +x_train_bool = Tensor.random(0_u8...2_u8, [bsz * 100, 2]) + +y_bool = x_train_bool[..., ...1] ^ x_train_bool[..., 1...] + +x_train = ctx.variable(x_train_bool.as_type(Float64)) +y = y_bool.as_type(Float64) + +net = Num::NN::Network.new(ctx) do + linear 2, 3 + relu + linear 3, 1 + sgd 0.7 + sigmoid_cross_entropy_loss +end + +losses = [] of Float64 + +50.times do |epoch| + 100.times do |batch_id| + offset = batch_id * 32 + x = x_train[offset...offset + 32] + target = y[offset...offset + 32] + + y_pred = net.forward(x) + + loss = net.loss(y_pred, target) + + puts "Epoch is: #{epoch}" + puts "Batch id: #{batch_id}" + puts "Loss is: #{loss.value.value}" + losses << loss.value.value + + loss.backprop + net.optimizer.update + end +end +``` + +``` +... +Epoch is: 49 +Batch id: 95 +Loss is: 0.00065050072686102 +Epoch is: 49 +Batch id: 96 +Loss is: 0.0006024037564266797 +Epoch is: 49 +Batch id: 97 +Loss is: 0.0005297538443899917 +Epoch is: 49 +Batch id: 98 +Loss is: 0.0005765025171222869 +Epoch is: 49 +Batch id: 99 +Loss is: 0.0005290653040218895 +``` + +The Network learns this function very quickly, as XOR is one of the simplest +distributions to hit. Since the training data is so limited, accuracy +can be a bit jagged, but eventually the network smooths out. + +### Loss over time + +![xorloss](xor_classifier_loss.png) diff --git a/examples/basic_xor_classifier/xor.cr b/examples/basic_xor_classifier/xor.cr new file mode 100644 index 00000000..668b41b5 --- /dev/null +++ b/examples/basic_xor_classifier/xor.cr @@ -0,0 +1,72 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +require "../../src/num" + +ctx = Num::Grad::Context(Tensor(Float64)).new + +bsz = 32 + +x_train_bool = Tensor.random(0_u8...2_u8, [bsz * 100, 2]) + +y_bool = x_train_bool[..., ...1] ^ x_train_bool[..., 1...] + +x_train = ctx.variable(x_train_bool.as_type(Float64)) +y = y_bool.as_type(Float64) + +net = Num::NN::Network.new(ctx) do + linear 2, 3 + relu + linear 3, 1 + sgd 0.7 + sigmoid_cross_entropy_loss +end + +losses = [] of Float64 + +50.times do |epoch| + 100.times do |batch_id| + offset = batch_id * 32 + x = x_train[offset...offset + 32] + target = y[offset...offset + 32] + + y_pred = net.forward(x) + + loss = net.loss(y_pred, target) + + puts "Epoch is: #{epoch}" + puts "Batch id: #{batch_id}" + puts "Loss is: #{loss.value.value}" + losses << loss.value.value + + loss.backprop + net.optimizer.update + end +end + +Num::Plot::Plot.plot do + scatter (0...losses.size), losses + x_label "Epochs" + y_label "Loss" + label "XOR Classifier Loss" +end diff --git a/examples/basic_xor_classifier/xor_classifier_loss.png b/examples/basic_xor_classifier/xor_classifier_loss.png new file mode 100644 index 00000000..abbd1ed5 Binary files /dev/null and b/examples/basic_xor_classifier/xor_classifier_loss.png differ diff --git a/examples/simple_scatter_plot/README.md b/examples/simple_scatter_plot/README.md new file mode 100644 index 00000000..bfc22d11 --- /dev/null +++ b/examples/simple_scatter_plot/README.md @@ -0,0 +1,15 @@ +## Simple Scatter Plot + +Using [PLplot](http://plplot.sourceforge.net/index.php), a fantastic library for scientific plotting, charts are extremely easy to create with Num.cr. + +```crystal +x = (0...100) +y = Tensor(Float64).rand([100]) + +Num::Plot::Plot.plot do + term nil + scatter x, y, code: 14, color: 2 +end +``` + +![scatter](simple_scatter.png) diff --git a/examples/simple_scatter_plot/scatter.cr b/examples/simple_scatter_plot/scatter.cr new file mode 100644 index 00000000..483903a2 --- /dev/null +++ b/examples/simple_scatter_plot/scatter.cr @@ -0,0 +1,32 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +require "../../src/num" + +x = (0...100) +y = Tensor(Float64).rand([100]) + +Num::Plot::Plot.plot do + term nil + scatter x, y, code: 14, color: 2 +end diff --git a/examples/simple_scatter_plot/simple_scatter.png b/examples/simple_scatter_plot/simple_scatter.png new file mode 100644 index 00000000..ded617e8 Binary files /dev/null and b/examples/simple_scatter_plot/simple_scatter.png differ diff --git a/shard.yml b/shard.yml index 6872f735..774df9f3 100644 --- a/shard.yml +++ b/shard.yml @@ -11,6 +11,8 @@ license: MIT dependencies: opencl: github: crystal-data/opencl.cr + alea: + github: nin93/alea scripts: postinstall: make ext diff --git a/spec/tensor/reduction_spec.cr b/spec/tensor/reduction_spec.cr new file mode 100644 index 00000000..3da33644 --- /dev/null +++ b/spec/tensor/reduction_spec.cr @@ -0,0 +1,31 @@ +require "../spec_helper" + +describe Tensor do + it "sorts a one dimensional Tensor" do + a = [4, 3, 2, 1].to_tensor + result = Num.sort(a) + expected = [1, 2, 3, 4].to_tensor + assert_array_equal(result, expected) + end + + it "sorts a strided Tensor" do + a = [4, 3, 2, 1].to_tensor[{..., 2}] + result = Num.sort(a) + expected = [2, 4] + assert_array_equal(result, expected) + end + + it "sorts a Tensor along an axis" do + a = [[3, 5, 6], [1, 1, 2], [9, 2, 3]].to_tensor + result = Num.sort(a, 0) + expected = [[1, 1, 2], [3, 2, 3], [9, 5, 6]].to_tensor + assert_array_equal(result, expected) + end + + it "sorts a strided Tensor along an axis" do + a = [[3, 4, 5, 1], [2, 1, 3, 2], [4, 7, 6, 2]].to_tensor[..., {..., 2}] + result = Num.sort(a, 0) + expected = [[2, 3], [3, 5], [4, 6]].to_tensor + assert_array_equal(result, expected) + end +end diff --git a/src/api.cr b/src/api.cr index 03e7b055..6a4ce468 100644 --- a/src/api.cr +++ b/src/api.cr @@ -1,5 +1,6 @@ require "./tensor/build" require "./tensor/creation" +require "./tensor/random" require "./tensor/linalg" require "./tensor/operators" require "./tensor/reductions" @@ -10,11 +11,27 @@ require "./cl_tensor/cl_tensor" require "./cl_tensor/creation" require "./cl_tensor/linalg" -require "./frame/frame" -require "./frame/series" -require "./frame/index" - require "./scikit/matrices" require "./scikit/clustering/kmeans" require "./libs/local" +require "./libs/nnpack" +require "./libs/plplot" + +require "./grad/primitives/*" +require "./grad/gates_arithmetic" +require "./grad/gates_blas" +require "./grad/variable_ops" + +require "./nn/primitives/*" +require "./nn/layers/*" +require "./nn/gates/*" +require "./nn/optimizer" +require "./nn/loss" +require "./nn/network" + +require "./nn/datasets/*" + +require "./plot/internal/*" +require "./plot/figures/*" +require "./plot/plot" diff --git a/src/frame/frame.cr b/src/frame/frame.cr deleted file mode 100644 index bb8a849d..00000000 --- a/src/frame/frame.cr +++ /dev/null @@ -1,256 +0,0 @@ -# Copyright (c) 2020 Crystal Data Contributors -# -# MIT License -# -# Permission is hereby granted, free of charge, to any person obtaining -# a copy of this software and associated documentation files (the -# "Software"), to deal in the Software without restriction, including -# without limitation the rights to use, copy, modify, merge, publish, -# distribute, sublicense, and/or sell copies of the Software, and to -# permit persons to whom the Software is furnished to do so, subject to -# the following conditions: -# -# The above copyright notice and this permission notice shall be -# included in all copies or substantial portions of the Software. -# -# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF -# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE -# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION -# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION -# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. - -require "../api" -require "./series" -require "./frame_slice" -require "csv" - -class DataFrame(T, V) - getter c : T - getter index : Index(V) - getter size : Int32 - - # Create a `DataFrame` from an index, and a NamedTuple - # of values. This is private so that all other methods - # can coerce the values of `data` to `Series`, so that - # no other types can be present in a `DataFrame` - # - # Arguments - # --------- - # `index` : Index(V) - # Index for all Series in the DataFrame - # `data` : NamedTuple - # Schema and Series of the DataFrame - def initialize(index : Index(V), size : Int, **data : **T) - @c = data - @index = index - @size = size.to_i - end - - # Create a `DataFrame` from a variadic number of arguments. - # - # The arguments can be of a flexible type, as long as they are - # one dimensional and can be cast to `Tensor`s, and therefore - # `Series`. - # - # Arguments - # --------- - # `data` : NamedTuple of Tensor's or Enumerables - # Data to convert into a `DataFrame` - # - # Examples - # -------- - # ``` - # a = Tensor.random(0.0...10.0, [4]) - # b = [1, 2, 3, 4] - # - # df = DataFrame.from_items(a: a, b: b) - # puts df - # - # a b - # 0 4.06169 1 - # 1 7.55353 2 - # 2 1.26119 3 - # 3 1.16003 4 - # ``` - def self.from_items(**data : **U) forall U - {% begin %} - data = NamedTuple.new( - {% for key, value in U %} - {{key}}: Series.new( - data[{{key.symbolize}}].to_tensor, - name: {{key.symbolize}} - ), - {% end %} - ) - {% end %} - size = self.check_attributes(data) - ii = Index.range(size) - data.each do |k, v| - v.set_index!(ii) - end - new(ii, size, **data) - end - - # :nodoc: - private def self.check_attributes(data) : Int32 - s0 = 0 - data.each_with_index do |k, v, i| - if i == 0 - s0 = v.size - end - - if v.size != s0 - raise Num::Internal::ShapeError.new("All inputs must be the same size") - end - end - s0 - end - - def [](i : V) - {% begin %} - data = NamedTuple.new( - {% for key, value in T %} - {{key}}: @c[{{key.symbolize}}][i], - {% end %} - ) - FrameSlice.new(**data) - {% end %} - end - - # :nodoc: - def each - @size.times do |i| - {% begin %} - data = NamedTuple.new( - {% for key, value in T %} - {{key}}: @c[{{key.symbolize}}].iat(i), - {% end %} - ) - yield data - {% end %} - end - end - - def each_with_index - i = 0 - each do |e| - yield e, @index.iat(i).key - i += 1 - end - end - - # :nodoc: - def to_s(io) - iw = @index.max_repr_width - dw = Hash(Symbol, Int32).new - @c.each do |k, v| - dw[k] = {"#{k}".size, v.max_repr_width}.max - end - io << " ".rjust(iw) - io << " " - @c.each do |k, v| - io << "#{k}".rjust(dw[k]) - io << " " - end - io << "\n" - i = 0 - each_with_index do |e, i| - io << "#{Num::Internal.format(i)}".rjust(iw) - io << " " - e.each do |k, v| - if k != :index - io << "#{Num::Internal.format(v)}".rjust(dw[k]) - io << " " - end - end - io << "\n" - end - end - - # :nodoc: - def to_csv - CSV.build do |csv| - csv.row do |r| - r << "index" - @c.each_key do |k| - r << k - end - end - each_with_index do |e, i| - csv.row do |r| - r << i - e.each_value do |v| - r << v - end - end - end - end - end - - # :nodoc: - macro reduce(reduction) - def {{reduction.id}} - \{% begin %} - FrameSlice.new( - \{% for key, value in T %} - \{{key}}: @c[\{{key.symbolize}}].{{reduction.id}}, - \{% end %} - ) - \{% end %} - end - end - - reduce sum - reduce prod - reduce min - reduce max - reduce mean - reduce std - reduce all - reduce any - reduce unique - reduce to_a - - # :nodoc: - macro elementwise(fn) - def {{fn.id}}(other) - \{% begin %} - data = NamedTuple.new( - \{% for key, value in T %} - \{{key}}: @c[\{{key.symbolize}}].{{fn.id}}(other), - \{% end %} - ) - DataFrame.new(@index, @size, **data) - \{% end %} - end - end - - elementwise add - elementwise subtract - elementwise multiply - elementwise divide - elementwise floordiv - elementwise power - elementwise modulo - elementwise left_shift - elementwise right_shift - elementwise bitwise_and - elementwise bitwise_or - elementwise bitwise_xor - elementwise equal - elementwise not_equal - elementwise greater - elementwise greater_equal - elementwise less - elementwise less_equal -end - -macro agg(nt, **fns) - FrameSlice.new( - {% for k, v in fns %} - {{k}}: {{nt}}.c[{{k.symbolize}}].{{v.id}}, - {% end %} - ) -end diff --git a/src/frame/index.cr b/src/frame/index.cr deleted file mode 100644 index b2ebba52..00000000 --- a/src/frame/index.cr +++ /dev/null @@ -1,107 +0,0 @@ -# Copyright (c) 2020 Crystal Data Contributors -# -# MIT License -# -# Permission is hereby granted, free of charge, to any person obtaining -# a copy of this software and associated documentation files (the -# "Software"), to deal in the Software without restriction, including -# without limitation the rights to use, copy, modify, merge, publish, -# distribute, sublicense, and/or sell copies of the Software, and to -# permit persons to whom the Software is furnished to do so, subject to -# the following conditions: -# -# The above copyright notice and this permission notice shall be -# included in all copies or substantial portions of the Software. -# -# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF -# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE -# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION -# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION -# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. - -# :nodoc: -class Hash(K, V) - def entries - @entries - end -end - -# An `Index` implements a `Hash` that allows rows to be quickly -# looked up in a `Frame` or `Series` in constant time. -# -# Many types of indexes will be supported, but currently only -# `Integer` indexes are implemented. -class Index(T) - @data : Hash(T, Int32) - - # Initialize an `Index` from a `Hash`. This should - # only ever be used by internal methods - private def initialize(@data : Hash(T, Int32)) - end - - # Create a range index from a size. This index will always - # monotonically increase and be unique. - # - # Arguments - # --------- - # `n` : Int - # Size of the index - # - # Examples - # -------- - # ``` - # Index.range(5) # => Index[0, 1, 2, 3, 4] - # ``` - def self.range(n : Int) - d = Hash(Int32, Int32).new - n.times do |i| - d[i] = i - end - new(d) - end - - # Find the corresponding row value of an index at a provided - # index value. - # - # Arguments - # --------- - # `i` : T - # Key to lookup - # - # Examples - # ``` - # i = Index.range(5) - # i[0] # => 0 - # ``` - def [](i : T) - @data[i] - end - - # Returns the appropriate `Entry` of an index at a given - # integer value. This relies on the fact that Crystal - # hashes are insertion sorted. - def iat(index : Int) - @data.entries[index] - end - - # :nodoc: - def each - @data.each do |k, v| - yield k, v - end - end - - # :nodoc: - def max_repr_width - w = 0 - each do |k, _| - l = Num::Internal.format(k).size - if l > w - w = l - end - end - w - end -end diff --git a/src/frame/series.cr b/src/frame/series.cr deleted file mode 100644 index 3e891bea..00000000 --- a/src/frame/series.cr +++ /dev/null @@ -1,202 +0,0 @@ -# Copyright (c) 2020 Crystal Data Contributors -# -# MIT License -# -# Permission is hereby granted, free of charge, to any person obtaining -# a copy of this software and associated documentation files (the -# "Software"), to deal in the Software without restriction, including -# without limitation the rights to use, copy, modify, merge, publish, -# distribute, sublicense, and/or sell copies of the Software, and to -# permit persons to whom the Software is furnished to do so, subject to -# the following conditions: -# -# The above copyright notice and this permission notice shall be -# included in all copies or substantial portions of the Software. -# -# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF -# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE -# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION -# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION -# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. - -require "./index" - -# A `Series` is a one dimensional view of a `Frame`. -# Slicing column-wise will return a `Series`. A `Series` -# stores its own index information even if an index -# is being tracked by a Frame. -class Series(T, V) - getter data : Tensor(T) - getter name : Symbol - getter index : Index(V) - - # Initializes a `Series` from its components. This - # method is only used by internal calls - # - # Arguments - # `data` : Tensor(T) - # Tensor containing the data of a series - # `name` : Symbol - # Series description - # `index` - # Hash allowing fast access to `Series` values - private def initialize(data : Tensor(T), name : Symbol, index : Index(V)) - @data = data - @name = name - @index = index - end - - # Initializes a `Series` from a `Tensor` and a `name`. - # If no name is provided the `Series` will be unnamed, - # and a name will be inferred if adding this to a - # `Frame` - # - # Arguments - # --------- - # `t` : Tensor(U) - # One dimensional input `Tensor` - # `name` : Symbol - # Identifier for the `Series` - def self.new(t : Tensor(U), name : Symbol = :unnamed) forall U - index = Index.range(t.size) - new(t, name, index) - end - - # Initializes a `Series` from a `Tensor`, `Index` and a `name`. - # If no name is provided the `Series` will be unnamed, - # and a name will be inferred if adding this to a - # `Frame` - # - # Arguments - # --------- - # `t` : Tensor(U) - # One dimensional input `Tensor` - # `index` : Index(V) - # Index for the `Series` - # `name` : Symbol - # Identifier for the `Series` - def self.new( - t : Tensor(U), - index : Index(V), - name : Symbol = :unnamed - ) forall U, V - new(t, name, index) - end - - # :nodoc: - def each - @data.each do |e| - yield e - end - end - - # :nodoc: - def each_with_index - @index.each do |k, v| - yield @data[v].value, k - end - end - - def [](i : V) - row = @index[i] - @data[row].value - end - - # :nodoc: - def iat(i : Int) - @data[i].value - end - - # :nodoc: - def size - @data.size - end - - def set_index!(index : Index(V)) - @index = index - end - - # :nodoc: - def to_s(io) - iw = @index.max_repr_width - vw = self.max_repr_width - each_with_index do |e, i| - io << "#{Num::Internal.format(i)}".ljust(iw) - io << " " - io << "#{Num::Internal.format(e)}".rjust(vw) - io << "\n" - end - io << "Name: #{@name}\n" - io << "dtype: #{T}" - end - - # :nodoc: - def max_repr_width - w = 0 - each do |k| - l = Num::Internal.format(k).size - if l > w - w = l - end - end - w - end - - private macro reduce(fn) - def {{fn.id}} - Num.{{fn.id}}(@data) - end - end - - reduce sum - reduce prod - reduce min - reduce max - reduce mean - reduce std - reduce all - reduce any - - private macro reduce_on(fn) - def {{fn.id}} - @data.{{fn.id}} - end - end - - reduce_on unique - reduce_on to_a - reduce_on to_tensor - - private macro elementwise(fn) - def {{fn.id}}(other : Number) - t = Num.{{fn.id}}(@data, other) - Series.new(t, @index, @name) - end - - def {{fn.id}}(other : FrameSlice) - t = Num.{{fn.id}}(@data, other.c[@name]) - Series.new(t, @index, @name) - end - end - - elementwise add - elementwise subtract - elementwise multiply - elementwise divide - elementwise floordiv - elementwise power - elementwise modulo - elementwise left_shift - elementwise right_shift - elementwise bitwise_and - elementwise bitwise_or - elementwise bitwise_xor - elementwise equal - elementwise not_equal - elementwise greater - elementwise greater_equal - elementwise less - elementwise less_equal -end diff --git a/src/grad/attribution.txt b/src/grad/attribution.txt new file mode 100644 index 00000000..c7a3b60d --- /dev/null +++ b/src/grad/attribution.txt @@ -0,0 +1,4 @@ +Most of the work in this module was inspired by @mratsim's fantastic +Arraymancer library. + +https://github.com/mratsim/Arraymancer diff --git a/src/grad/gates_arithmetic.cr b/src/grad/gates_arithmetic.cr new file mode 100644 index 00000000..df6e63a8 --- /dev/null +++ b/src/grad/gates_arithmetic.cr @@ -0,0 +1,145 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# :nodoc: +class Num::Grad::AddGate(T) < Num::Grad::Gate(T) + # :nodoc: + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [gradient, gradient] + end + + # :nodoc: + def cache(result : Num::Grad::Variable(T), *args) + a, b = args + + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Add", self, result, a, b) + end +end + +# :nodoc: +class Num::Grad::SubtractGate(T) < Num::Grad::Gate(T) + # :nodoc: + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [gradient, -gradient] + end + + # :nodoc: + def cache(result : Num::Grad::Variable(T), *args) + a, b = args + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Sub", self, result, a, b) + end +end + +# :nodoc: +class Num::Grad::MultiplyGate(T) < Num::Grad::Gate(T) + getter a : Num::Grad::Variable(T) + getter b : Num::Grad::Variable(T) + + # :nodoc: + def initialize(@a : Num::Grad::Variable(T), @b : Num::Grad::Variable(T)) + end + + # :nodoc: + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + [gradient * @b.value, @a.value * gradient] + end + + # :nodoc: + def cache(result : Num::Grad::Variable(T), *args) + a, b = args + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Mul", self, result, a, b) + end +end + +# :nodoc: +class Num::Grad::DivideGate(T) < Num::Grad::Gate(T) + getter a : Num::Grad::Variable(T) + getter b : Num::Grad::Variable(T) + + # :nodoc: + def initialize(@a : Num::Grad::Variable(T), @b : Num::Grad::Variable(T)) + end + + # :nodoc: + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + r0 = gradient.map(@b.value) { |i, j| i / j } + r1 = gradient.map(@a.value, @b.value) { |i, j, k| -i * j / (k ** 2) } + [r0, r1] + end + + # :nodoc: + def cache(result : Num::Grad::Variable(T), *args) + a, b = args + result.grad = T.zeros_like(result.value) + result.requires_grad = true + Num::Grad.register("Div", self, result, a, b) + end +end + +# :nodoc: +class Num::Grad::PowerGate(T) < Num::Grad::Gate(T) + getter a : Num::Grad::Variable(T) + getter b : Num::Grad::Variable(T) + + # :nodoc: + def initialize(@a : Num::Grad::Variable(T), @b : Num::Grad::Variable(T)) + end + + # :nodoc: + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + r0 = gradient.map(a.value, b.value) do |grad, x, y| + grad * y * (x ** (y == 0 ? 1 : y - 1)) + end + + r1 = gradient.map(a.value, b.value) do |grad, x, y| + grad * (x ** y) * Math.log(x == 0 ? 1 : x) + end + + [r0, r1] + end + + # :nodoc: + def cache(result : Num::Grad::Variable(T), *args) + a, b = args + result.grad = T.zeros_like(result.value) + result.requires_grad = true + Num::Grad.register("Pow", self, result, a, b) + end +end diff --git a/src/grad/gates_blas.cr b/src/grad/gates_blas.cr new file mode 100644 index 00000000..3e13c711 --- /dev/null +++ b/src/grad/gates_blas.cr @@ -0,0 +1,52 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# :nodoc: +class Num::Grad::MatMulGate(T) < Num::Grad::Gate(T) + getter a : Num::Grad::Variable(T) + getter b : Num::Grad::Variable(T) + + # :nodoc: + def initialize(@a : Num::Grad::Variable(T), @b : Num::Grad::Variable(T)) + end + + # :nodoc: + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + r0 = gradient.matmul(@b.value.transpose) + r1 = @a.value.transpose.matmul(gradient) + + [r0, r1] + end + + # :nodoc: + def cache(result : Num::Grad::Variable(T), *args : Num::Grad::Variable(T)) + a, b = args + + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("MatMul", self, result, a, b) + end +end diff --git a/src/grad/primitives/context.cr b/src/grad/primitives/context.cr new file mode 100644 index 00000000..4425ba6d --- /dev/null +++ b/src/grad/primitives/context.cr @@ -0,0 +1,147 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# A Context keeps track of the computational graph for +# a number of operations. Variables that interact with each +# other must belong to the same context, or state will be +# lost while tracking operations done. +# +# The generic type of a context is always going to be a specific +# class of `Tensor`, to allow easy creation of gradients on the +# fly. Unlike standard `Tensor` operations, a `Context` cannot +# shift it's generic type, and operations resulting in a different +# data type will raise. + +class Num::Grad::Context(T) + # A list of all variables present in an operation. + # This list can contain duplicates + getter nodes : Array(Num::Grad::Node(T)) + + # If no_grad is set to true, operations will not + # be cached, and backpropogation will not be possible + getter no_grad : Bool + + # Contexts can only be initialized as empty, and + # a generic type must be provided + def initialize + @nodes = Array(Num::Grad::Node(T)).new + @no_grad = false + end + + # :nodoc: + def size : Int32 + @nodes.size + end + + # :nodoc: + def <<(node : Num::Grad::Node(T)) + @nodes << node + end + + def last : Num::Grad::Node(T) + # :nodoc: + @nodes.last + end + + # :nodoc: + def pop : Num::Grad::Node(T) + @nodes.pop + end + + # Creates a new variable within the `Context`. This variable + # must be able to be cast to a `Tensor` of type `T`. + # + # Arguments + # --------- + # *value* : Tensor-like + # A value that can be converted to a `Tensor` + # *requires_grad* : Bool + # Flag to indicate if operations should be cached for this + # variable + # + # Examples + # -------- + # ``` + # ctx = Context(Tensor(Float64)).new + # ctx.variable([1.0, 2.0, 3.0]) + # ``` + def variable(value, requires_grad : Bool = true) : Num::Grad::Variable(T) + Num::Grad::Variable.new(self, value.to_tensor, requires_grad) + end + + # :nodoc: + def to_s(io) + @nodes.each_with_index do |node, i| + io << "Node #{i}: #{node.name} - " + if node.parents.size <= 1 + io << node.parents[0].value.shape + else + io << "(" + node.parents.each_with_index do |parent, pi| + if pi != 0 + io << ", " + end + io << parent.value.shape + end + io << ")" + end + io << node.payload.variable.value.shape + unless i == @nodes.size - 1 + io << "\n" + end + end + end +end + +module Num::Grad + extend self + + # Cached a node in the computational graph. This is only required + # if an operation needs to be backpropogated. + # + # Arguments + # --------- + # *name* : String + # Description of the operation + # *gate* : Gate(U) + # Operation gate containing a backward method and + # cached arguments + # *result* : Variable(U) + # The result of the operation being cached + # *parents* : *Variable(U) + # The operands present in operation being cached + # + # This method should be used sparingly by Users of the application. + # It should only be necessary when a User is defining their own + # custom activation function or Layer. + def register( + name : String, + gate : Num::Grad::Gate(U), + result : Num::Grad::Variable(U), + *parents : Num::Grad::Variable(U) + ) forall U + payload = Num::Grad::Payload.new(result) + node = Num::Grad::Node.new(gate, parents.to_a, payload, name) + parents[0].context << node + end +end diff --git a/src/grad/primitives/gate.cr b/src/grad/primitives/gate.cr new file mode 100644 index 00000000..0066a547 --- /dev/null +++ b/src/grad/primitives/gate.cr @@ -0,0 +1,38 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# A Gate is an object that can cache the result of an operation, +# as well as backpropogate a payload backwards along the +# computational graph +# +# Child classes that inherit from this class can add instance +# variables if additional caching is needed, and these need +# to be populated when writing the cached operation +abstract class Num::Grad::Gate(T) + # Propogates an operation backwards, transforming a payload + # and returning an array of Tensors + abstract def backward(payload : Num::Grad::Payload(T)) : Array(T) + + # Caches the result of an operation on a context + abstract def cache(result : Num::Grad::Variable(T), *args) +end diff --git a/src/grad/primitives/node.cr b/src/grad/primitives/node.cr new file mode 100644 index 00000000..91696449 --- /dev/null +++ b/src/grad/primitives/node.cr @@ -0,0 +1,53 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# A Node is a member of a computational graph that contains +# a reference to a gate, as well as the parents of the operation +# and the payload that resulted from the operation. +class Num::Grad::Node(T) + # A Gate containing a backwards and cache function for + # a node + getter gate : Num::Grad::Gate(T) + + # The variables that created this node + getter parents : Array(Num::Grad::Variable(T)) + + # Wrapper around a Tensor, contains operation data + getter payload : Num::Grad::Payload(T) + + # Debug use only, contains a name for a node + getter name : String + + # This initializer shouldn't ever be called outside of + # Num::Grad.register. + # + # Users defining custom gradients and gates should + # follow the same rule + def initialize( + @gate : Gate(T), + @parents : Array(Num::Grad::Variable(T)), + @payload : Num::Grad::Payload(T), + @name : String = "" + ) + end +end diff --git a/src/grad/primitives/payload.cr b/src/grad/primitives/payload.cr new file mode 100644 index 00000000..af063749 --- /dev/null +++ b/src/grad/primitives/payload.cr @@ -0,0 +1,34 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# A Payload is a simple wrapper around a Variable. It +# is only abstracted out to be a bit more explicit that +# it is being passed around through an operation +class Num::Grad::Payload(T) + # Contents of the paylod + getter variable : Num::Grad::Variable(T) + + # This should only be called by internal methods + def initialize(@variable : Num::Grad::Variable(T)) + end +end diff --git a/src/grad/primitives/variable.cr b/src/grad/primitives/variable.cr new file mode 100644 index 00000000..a66ccd79 --- /dev/null +++ b/src/grad/primitives/variable.cr @@ -0,0 +1,107 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# A variable is an abstraction of a Tensor that tracks +# the operations done to the Tensor. It also keeps +# track of the gradient of the operation if a Variable +# needs to backpropogate. +# +# This is the fundamental object used in automatic +# differentiation, as well as the neural network aspects +# of Num.cr +class Num::Grad::Variable(T) + # The graph the variable is associated with. This is a reference, + # as a variable does not own its context + getter context : Num::Grad::Context(T) + + # The value of the Variable. This should not be edited outside + # of Variable operations, as other edits will not be tracked + # and will lead to incorrect results + getter value : T + + # The gradient of the Variable. This is set as a reference to + # the value of a Variable unless `backprop` has been called, in + # which case all related Variables will have their gradient + # updated correctly + property grad : T + + # If set to true, this variable will track its operations, + # otherwise it will act similar to a Tensor, only calculating + # forward operations + property requires_grad : Bool + + # Initialization method for a Variable. + # + # This method should only be called by a context, as it creates + # a Variable. Context provides a helper method to add a + # Variable to the computational graph that handles ownership + # of the context and other related instance variables + def initialize( + @context : Num::Grad::Context(T), + @value : T, + @requires_grad : Bool = false + ) + if @requires_grad + @grad = T.zeros_like(@value) + else + @grad = @value + end + end + + # :nodoc: + def is_grad_needed : Bool + @requires_grad && !@context.no_grad + end + + # :nodoc: + def to_s(io) + @value.to_s(io) + end + + # Back propogates an operation along a computational graph. + # This operation will destroy the operational graph, populating + # the gradients for all variables that are predecessors of + # the Variable this is called on. + # + # Even if this is called on the first node in a graph, it will + # destroy all descendents of this variable stored by the + # Context + def backprop + @grad = T.ones_like(@value) + + while @context.size > 0 && @context.last.payload.variable != self + @context.pop + end + + while @context.size > 0 + cur_node = @context.pop + diffs = cur_node.gate.backward(cur_node.payload) + diffs.each_with_index do |diff, i| + parent_i = cur_node.parents[i] + if parent_i.requires_grad + parent_i.grad += diff + end + end + end + end +end diff --git a/src/grad/variable_ops.cr b/src/grad/variable_ops.cr new file mode 100644 index 00000000..ff1a6983 --- /dev/null +++ b/src/grad/variable_ops.cr @@ -0,0 +1,214 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# A variable is an abstraction of a Tensor that tracks +# the operations done to the Tensor. It also keeps +# track of the gradient of the operation if a Variable +# needs to backpropogate. +# +# This is the fundamental object used in automatic +# differentiation, as well as the neural network aspects +# of Num.cr +class Num::Grad::Variable(T) + private macro operator_op(operator, gate_cls, *args) + def {{operator.id}}(other : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + result = @context.variable(@value {{operator.id}} other.value) + + if self.is_grad_needed || other.is_grad_needed + gate = {{gate_cls}}.new {{*args}} + gate.cache(result, self, other) + end + + result + end + end + + # Adds a variable to another variable and stores + # the derivative of the operation in the computational + # graph. + # + # Arguments + # --------- + # *other* : Num::Grad::Variable(T) + # - right hand side of the operation + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([2.0]) + # b = ctx.variable([3.0]) + # + # f = a + b # => [5.0] + # f.backprop + # ``` + operator_op :+, Num::Grad::AddGate(T) + + # Subtracts a variable from another variable and stores + # the derivative of the operation in the computational + # graph. + # + # Arguments + # --------- + # *other* : Num::Grad::Variable(T) + # - right hand side of the operation + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([2.0]) + # b = ctx.variable([3.0]) + # + # f = a - b # => [-1.0] + # f.backprop + # ``` + operator_op :-, Num::Grad::SubtractGate(T) + + # Multiples a variable to another variable and stores + # the derivative of the operation in the computational + # graph. + # + # Arguments + # --------- + # *other* : Num::Grad::Variable(T) + # - right hand side of the operation + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([2.0]) + # b = ctx.variable([3.0]) + # + # f = a * b # => [6.0] + # f.backprop + # ``` + operator_op :*, Num::Grad::MultiplyGate(T), self, other + + # Raises a variable to another variable and stores + # the derivative of the operation in the computational + # graph. + # + # Arguments + # --------- + # *other* : Num::Grad::Variable(T) + # - right hand side of the operation + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([2.0]) + # b = ctx.variable([3.0]) + # + # f = a ** b # => [8.0] + # f.backprop + # ``` + operator_op :**, Num::Grad::PowerGate(T), self, other + + # Divides a variable by another variable and stores + # the derivative of the operation in the computational + # graph. + # + # Arguments + # --------- + # *other* : Num::Grad::Variable(T) + # - right hand side of the operation + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([2.0]) + # b = ctx.variable([3.0]) + # + # f = a / b # => [0.66667] + # f.backprop + # ``` + operator_op :/, Num::Grad::DivideGate(T), self, other + + # Matrix multiply operator for two variables. Computes the + # dot product of two matrices and stores the result in the + # computational graph + # + # Arguments + # --------- + # *other* : Num::Grad::Variable(T) + # - right hand side of the operation + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([[2.0], [2.0]]) + # b = ctx.variable([[3.0, 3.0]]) + # + # f = a + b + # + # # [[6, 6], + # # [6, 6]] + # + # f.backprop + # ``` + def matmul(b : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + result = Num::Grad::Variable.new(@context, self.value.matmul(b.value)) + + if self.is_grad_needed || b.is_grad_needed + gate = Num::Grad::MatMulGate(T).new(self, b) + gate.cache(result, self, b) + end + result + end + + # Slices a variable. Slices the gradient of the variable + # using the same arguments + # + # Arguments + # --------- + # *args* + # Slicing arguments, slicing behavior is the same as + # it is for a standard Tensor + # + # Examples + # -------- + # ``` + # ctx = Num::Grad::Context(Tensor(Float64)).new + # + # a = ctx.variable([[2.0], [3.0]]) + # b = a[1] + # b # => [3] + # ``` + def [](*args) + output = @value[*args] + result = @context.variable(output, requires_grad: @requires_grad) + result.grad = @grad[*args] + result + end +end diff --git a/src/libs/nnpack.cr b/src/libs/nnpack.cr new file mode 100644 index 00000000..949839a3 --- /dev/null +++ b/src/libs/nnpack.cr @@ -0,0 +1,162 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +@[Link("nnpack")] +lib LibNNPACK + enum NNPStatus : LibC::Int + SUCCESS = 0 + INVALID_BATCH_SIZE = 2 + INVALID_CHANNELS = 3 + INVALID_INPUT_CHANNELS = 4 + INVALID_OUTPUT_CHANNELS = 5 + INVALID_INPUT_SIZE = 10 + INVALID_INPUT_STRIDE = 11 + INVALID_INPUT_PADDING = 12 + INVALID_KERNEL_SIZE = 13 + INVALID_POOLING_SIZE = 14 + INVALID_POOLING_STRIDE = 15 + INVALID_ALGORITHM = 16 + INVALID_TRANSFORM_STRATEGY = 17 + UNSUPPORTED_INPUT_SIZE = 20 + UNSUPPORTED_INPUT_STRIDE = 21 + UNSUPPORTED_INPUT_PADDING = 22 + UNSUPPORTED_KERNEL_SIZE = 23 + UNSUPPORTED_POOLING_SIZE = 24 + UNSUPPORTED_POOLING_STRIDE = 25 + UNSUPPORTED_ALGORITHM = 26 + UNSUPPORTED_TRANSFORM_STRATEGY = 27 + UNSUPPORTED_ACTIVATION = 28 + UNSUPPORTED_ACTIVATION_PARAMETERS = 29 + UNINITIALIZED = 50 + UNSUPPORTED_HARDWARE = 51 + OUT_OF_MEMORY = 52 + INSUFFICIENT_BUFFER = 53 + MISALIGNED_BUFFER = 54 + end + + enum NNPActivation : LibC::Int + IDENTITY = 0 + RELU = 1 + end + + enum NNPConvolutionAlgorithm : LibC::Int + AUTO = 0 + FT8X8 = 1 + FT16X16 = 2 + WT8X8 = 3 + IMPLICIT_GEMM = 4 + DIRECT = 5 + WT8X8_FP16 = 6 + end + + enum NNPConvolutionTransformStrategy : LibC::Int + COMPUTE = 1 + PRECOMPUTE = 2 + REUSE = 3 + end + + struct NNPSize + width : LibC::SizeT + height : LibC::SizeT + end + + struct NNPPadding + top : LibC::SizeT + right : LibC::SizeT + bottom : LibC::SizeT + left : LibC::SizeT + end + + struct NNPProfile + total : LibC::Double + input_transform : LibC::Double + kernel_transform : LibC::Double + output_transform : LibC::Double + block_multiplication : LibC::Double + end + + fun initialize = nnp_initialize : NNPStatus + fun deinitialize = nnp_deinitialize : NNPStatus + + fun nnp_convolution_output = nnp_convolution_output( + algorithm : NNPConvolutionAlgorithm, + batch_size : LibC::SizeT, + input_channels : LibC::SizeT, + output_channels : LibC::SizeT, + input_size : NNPSize, + input_padding : NNPPadding, + kernel_size : NNPSize, + input : LibC::Float*, + kernel : LibC::Float*, + bias : LibC::Float*, + output : LibC::Float*, + workspace_buffer : Void*, + workspace_size : LibC::SizeT*, + activation : NNPActivation, + activation_parameters : Void*, + threadpool : Void*, + profile : NNPProfile* + ) : NNPStatus + + fun nnp_convolution_input_gradient = nnp_convolution_input_gradient( + algorithm : NNPConvolutionAlgorithm, + batch_size : LibC::SizeT, + input_channels : LibC::SizeT, + output_channels : LibC::SizeT, + input_size : NNPSize, + input_padding : NNPPadding, + kernel_size : NNPSize, + grad_output : LibC::Float*, + kernel : LibC::Float*, + grad_input : LibC::Float*, + workspace_buffer : Void*, + workspace_size : LibC::SizeT*, + activation : NNPActivation, + activation_parameters : Void*, + threadpool : Void*, + profile : NNPProfile* + ) : NNPStatus + + fun nnp_convolution_kernel_gradient = nnp_convolution_kernel_gradient( + algorithm : NNPConvolutionAlgorithm, + batch_size : LibC::SizeT, + input_channels : LibC::SizeT, + output_channels : LibC::SizeT, + input_size : NNPSize, + input_padding : NNPPadding, + kernel_size : NNPSize, + input : LibC::Float*, + grad_output : LibC::Float*, + grad_kernel : LibC::Float*, + workspace_buffer : Void*, + workspace_size : LibC::SizeT*, + activation : NNPActivation, + activation_parameters : Void*, + threadpool : Void*, + profile : NNPProfile* + ) : NNPStatus +end + +{% if flag?(:nnpack) %} + LibNNPACK.initialize +{% end %} diff --git a/src/libs/plplot.cr b/src/libs/plplot.cr new file mode 100644 index 00000000..8daad1d2 --- /dev/null +++ b/src/libs/plplot.cr @@ -0,0 +1,228 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +@[Link("plplot")] +lib LibPlplot + fun pl_setcontlabelformat = c_pl_setcontlabelformat(lexp : Plint, sigdig : Plint) + alias X__Int32T = LibC::Int + alias Int32T = X__Int32T + alias Plint = Int32T + fun pl_setcontlabelparam = c_pl_setcontlabelparam(offset : Plflt, size : Plflt, spacing : Plflt, active : Plint) + alias Plflt = LibC::Double + fun pladv = c_pladv(page : Plint) + fun plarc = c_plarc(x : Plflt, y : Plflt, a : Plflt, b : Plflt, angle1 : Plflt, angle2 : Plflt, rotate : Plflt, fill : Plbool) + alias Plbool = Plint + fun plaxes = c_plaxes(x0 : Plflt, y0 : Plflt, xopt : PlcharVector, xtick : Plflt, nxsub : Plint, yopt : PlcharVector, ytick : Plflt, nysub : Plint) + alias PlcharVector = LibC::Char* + fun plbin = c_plbin(nbin : Plint, x : PlfltVector, y : PlfltVector, opt : Plint) + alias PlfltVector = Plflt* + fun plbtime = c_plbtime(year : PlintNcScalar, month : PlintNcScalar, day : PlintNcScalar, hour : PlintNcScalar, min : PlintNcScalar, sec : PlfltNcScalar, ctime : Plflt) + alias PlintNcScalar = Plint* + alias PlfltNcScalar = Plflt* + fun plbop = c_plbop + fun plbox = c_plbox(xopt : PlcharVector, xtick : Plflt, nxsub : Plint, yopt : PlcharVector, ytick : Plflt, nysub : Plint) + fun plbox3 = c_plbox3(xopt : PlcharVector, xlabel : PlcharVector, xtick : Plflt, nxsub : Plint, yopt : PlcharVector, ylabel : PlcharVector, ytick : Plflt, nysub : Plint, zopt : PlcharVector, zlabel : PlcharVector, ztick : Plflt, nzsub : Plint) + fun plcalc_world = c_plcalc_world(rx : Plflt, ry : Plflt, wx : PlfltNcScalar, wy : PlfltNcScalar, window : PlintNcScalar) + fun plclear = c_plclear + fun plcol0 = c_plcol0(icol0 : Plint) + fun plcol1 = c_plcol1(col1 : Plflt) + fun plconfigtime = c_plconfigtime(scale : Plflt, offset1 : Plflt, offset2 : Plflt, ccontrol : Plint, ifbtime_offset : Plbool, year : Plint, month : Plint, day : Plint, hour : Plint, min : Plint, sec : Plflt) + fun plcont = c_plcont(f : PlfltMatrix, nx : Plint, ny : Plint, kx : Plint, lx : Plint, ky : Plint, ly : Plint, clevel : PlfltVector, nlevel : Plint, pltr : PltransformCallback, pltr_data : PlPointer) + alias PlfltMatrix = Plflt** + alias PlPointer = Void* + alias PltransformCallback = (Plflt, Plflt, PlfltNcScalar, PlfltNcScalar, PlPointer -> Void) + fun plcpstrm = c_plcpstrm(iplsr : Plint, flags : Plbool) + fun plctime = c_plctime(year : Plint, month : Plint, day : Plint, hour : Plint, min : Plint, sec : Plflt, ctime : PlfltNcScalar) + fun plend = c_plend + fun plend1 = c_plend1 + fun plenv = c_plenv(xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, just : Plint, axis : Plint) + fun plenv0 = c_plenv0(xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, just : Plint, axis : Plint) + fun pleop = c_pleop + fun plerrx = c_plerrx(n : Plint, xmin : PlfltVector, xmax : PlfltVector, y : PlfltVector) + fun plerry = c_plerry(n : Plint, x : PlfltVector, ymin : PlfltVector, ymax : PlfltVector) + fun plfamadv = c_plfamadv + fun plfill = c_plfill(n : Plint, x : PlfltVector, y : PlfltVector) + fun plfill3 = c_plfill3(n : Plint, x : PlfltVector, y : PlfltVector, z : PlfltVector) + fun plflush = c_plflush + fun plfont = c_plfont(ifont : Plint) + fun plfontld = c_plfontld(fnt : Plint) + fun plgchr = c_plgchr(p_def : PlfltNcScalar, p_ht : PlfltNcScalar) + fun plgcmap1_range = c_plgcmap1_range(min_color : PlfltNcScalar, max_color : PlfltNcScalar) + fun plgcol0 = c_plgcol0(icol0 : Plint, r : PlintNcScalar, g : PlintNcScalar, b : PlintNcScalar) + fun plgcol0a = c_plgcol0a(icol0 : Plint, r : PlintNcScalar, g : PlintNcScalar, b : PlintNcScalar, alpha : PlfltNcScalar) + fun plgcolbg = c_plgcolbg(r : PlintNcScalar, g : PlintNcScalar, b : PlintNcScalar) + fun plgcolbga = c_plgcolbga(r : PlintNcScalar, g : PlintNcScalar, b : PlintNcScalar, alpha : PlfltNcScalar) + fun plgcompression = c_plgcompression(compression : PlintNcScalar) + fun plgdev = c_plgdev(p_dev : PlcharNcVector) + alias PlcharNcVector = LibC::Char* + fun plgdidev = c_plgdidev(p_mar : PlfltNcScalar, p_aspect : PlfltNcScalar, p_jx : PlfltNcScalar, p_jy : PlfltNcScalar) + fun plgdiori = c_plgdiori(p_rot : PlfltNcScalar) + fun plgdiplt = c_plgdiplt(p_xmin : PlfltNcScalar, p_ymin : PlfltNcScalar, p_xmax : PlfltNcScalar, p_ymax : PlfltNcScalar) + fun plgdrawmode = c_plgdrawmode : Plint + fun plgfci = c_plgfci(p_fci : PlunicodeNcScalar) + alias X__Uint32T = LibC::UInt + alias Uint32T = X__Uint32T + alias Pluint = Uint32T + alias Plunicode = Pluint + alias PlunicodeNcScalar = Plunicode* + fun plgfam = c_plgfam(p_fam : PlintNcScalar, p_num : PlintNcScalar, p_bmax : PlintNcScalar) + fun plgfnam = c_plgfnam(fnam : PlcharNcVector) + fun plgfont = c_plgfont(p_family : PlintNcScalar, p_style : PlintNcScalar, p_weight : PlintNcScalar) + fun plglevel = c_plglevel(p_level : PlintNcScalar) + fun plgpage = c_plgpage(p_xp : PlfltNcScalar, p_yp : PlfltNcScalar, p_xleng : PlintNcScalar, p_yleng : PlintNcScalar, p_xoff : PlintNcScalar, p_yoff : PlintNcScalar) + fun plgra = c_plgra + fun plgradient = c_plgradient(n : Plint, x : PlfltVector, y : PlfltVector, angle : Plflt) + fun plgriddata = c_plgriddata(x : PlfltVector, y : PlfltVector, z : PlfltVector, npts : Plint, xg : PlfltVector, nptsx : Plint, yg : PlfltVector, nptsy : Plint, zg : PlfltNcMatrix, type : Plint, data : Plflt) + alias PlfltNcMatrix = Plflt** + fun plgspa = c_plgspa(xmin : PlfltNcScalar, xmax : PlfltNcScalar, ymin : PlfltNcScalar, ymax : PlfltNcScalar) + fun plgstrm = c_plgstrm(p_strm : PlintNcScalar) + fun plgver = c_plgver(p_ver : PlcharNcVector) + fun plgvpd = c_plgvpd(p_xmin : PlfltNcScalar, p_xmax : PlfltNcScalar, p_ymin : PlfltNcScalar, p_ymax : PlfltNcScalar) + fun plgvpw = c_plgvpw(p_xmin : PlfltNcScalar, p_xmax : PlfltNcScalar, p_ymin : PlfltNcScalar, p_ymax : PlfltNcScalar) + fun plgxax = c_plgxax(p_digmax : PlintNcScalar, p_digits : PlintNcScalar) + fun plgyax = c_plgyax(p_digmax : PlintNcScalar, p_digits : PlintNcScalar) + fun plgzax = c_plgzax(p_digmax : PlintNcScalar, p_digits : PlintNcScalar) + fun plhist = c_plhist(n : Plint, data : PlfltVector, datmin : Plflt, datmax : Plflt, nbin : Plint, opt : Plint) + fun plhlsrgb = c_plhlsrgb(h : Plflt, l : Plflt, s : Plflt, p_r : PlfltNcScalar, p_g : PlfltNcScalar, p_b : PlfltNcScalar) + fun plinit = c_plinit + fun pljoin = c_pljoin(x1 : Plflt, y1 : Plflt, x2 : Plflt, y2 : Plflt) + fun pllab = c_pllab(xlabel : PlcharVector, ylabel : PlcharVector, tlabel : PlcharVector) + fun pllegend = c_pllegend(p_legend_width : PlfltNcScalar, p_legend_height : PlfltNcScalar, opt : Plint, position : Plint, x : Plflt, y : Plflt, plot_width : Plflt, bg_color : Plint, bb_color : Plint, bb_style : Plint, nrow : Plint, ncolumn : Plint, nlegend : Plint, opt_array : PlintVector, text_offset : Plflt, text_scale : Plflt, text_spacing : Plflt, text_justification : Plflt, text_colors : PlintVector, text : PlcharMatrix, box_colors : PlintVector, box_patterns : PlintVector, box_scales : PlfltVector, box_line_widths : PlfltVector, line_colors : PlintVector, line_styles : PlintVector, line_widths : PlfltVector, symbol_colors : PlintVector, symbol_scales : PlfltVector, symbol_numbers : PlintVector, symbols : PlcharMatrix) + alias PlintVector = Plint* + alias PlcharMatrix = LibC::Char** + fun plcolorbar = c_plcolorbar(p_colorbar_width : PlfltNcScalar, p_colorbar_height : PlfltNcScalar, opt : Plint, position : Plint, x : Plflt, y : Plflt, x_length : Plflt, y_length : Plflt, bg_color : Plint, bb_color : Plint, bb_style : Plint, low_cap_color : Plflt, high_cap_color : Plflt, cont_color : Plint, cont_width : Plflt, n_labels : Plint, label_opts : PlintVector, labels : PlcharMatrix, n_axes : Plint, axis_opts : PlcharMatrix, ticks : PlfltVector, sub_ticks : PlintVector, n_values : PlintVector, values : PlfltMatrix) + fun pllightsource = c_pllightsource(x : Plflt, y : Plflt, z : Plflt) + fun plline = c_plline(n : Plint, x : PlfltVector, y : PlfltVector) + fun plline3 = c_plline3(n : Plint, x : PlfltVector, y : PlfltVector, z : PlfltVector) + fun pllsty = c_pllsty(lin : Plint) + fun plmap = c_plmap(mapform : PlmapformCallback, name : PlcharVector, minx : Plflt, maxx : Plflt, miny : Plflt, maxy : Plflt) + alias PlfltNcVector = Plflt* + alias PlmapformCallback = (Plint, PlfltNcVector, PlfltNcVector -> Void) + fun plmapline = c_plmapline(mapform : PlmapformCallback, name : PlcharVector, minx : Plflt, maxx : Plflt, miny : Plflt, maxy : Plflt, plotentries : PlintVector, nplotentries : Plint) + fun plmapstring = c_plmapstring(mapform : PlmapformCallback, name : PlcharVector, string : PlcharVector, minx : Plflt, maxx : Plflt, miny : Plflt, maxy : Plflt, plotentries : PlintVector, nplotentries : Plint) + fun plmaptex = c_plmaptex(mapform : PlmapformCallback, name : PlcharVector, dx : Plflt, dy : Plflt, just : Plflt, text : PlcharVector, minx : Plflt, maxx : Plflt, miny : Plflt, maxy : Plflt, plotentry : Plint) + fun plmapfill = c_plmapfill(mapform : PlmapformCallback, name : PlcharVector, minx : Plflt, maxx : Plflt, miny : Plflt, maxy : Plflt, plotentries : PlintVector, nplotentries : Plint) + fun plmeridians = c_plmeridians(mapform : PlmapformCallback, dlong : Plflt, dlat : Plflt, minlong : Plflt, maxlong : Plflt, minlat : Plflt, maxlat : Plflt) + fun plmesh = c_plmesh(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint) + fun plmeshc = c_plmeshc(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint, clevel : PlfltVector, nlevel : Plint) + fun plmkstrm = c_plmkstrm(p_strm : PlintNcScalar) + fun plmtex = c_plmtex(side : PlcharVector, disp : Plflt, pos : Plflt, just : Plflt, text : PlcharVector) + fun plmtex3 = c_plmtex3(side : PlcharVector, disp : Plflt, pos : Plflt, just : Plflt, text : PlcharVector) + fun plot3d = c_plot3d(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint, side : Plbool) + fun plot3dc = c_plot3dc(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint, clevel : PlfltVector, nlevel : Plint) + fun plot3dcl = c_plot3dcl(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint, clevel : PlfltVector, nlevel : Plint, indexxmin : Plint, indexxmax : Plint, indexymin : PlintVector, indexymax : PlintVector) + fun plpat = c_plpat(nlin : Plint, inc : PlintVector, del : PlintVector) + fun plpath = c_plpath(n : Plint, x1 : Plflt, y1 : Plflt, x2 : Plflt, y2 : Plflt) + fun plpoin = c_plpoin(n : Plint, x : PlfltVector, y : PlfltVector, code : Plint) + fun plpoin3 = c_plpoin3(n : Plint, x : PlfltVector, y : PlfltVector, z : PlfltVector, code : Plint) + fun plpoly3 = c_plpoly3(n : Plint, x : PlfltVector, y : PlfltVector, z : PlfltVector, draw : PlboolVector, ifcc : Plbool) + alias PlboolVector = Plbool* + fun plprec = c_plprec(setp : Plint, prec : Plint) + fun plpsty = c_plpsty(patt : Plint) + fun plptex = c_plptex(x : Plflt, y : Plflt, dx : Plflt, dy : Plflt, just : Plflt, text : PlcharVector) + fun plptex3 = c_plptex3(wx : Plflt, wy : Plflt, wz : Plflt, dx : Plflt, dy : Plflt, dz : Plflt, sx : Plflt, sy : Plflt, sz : Plflt, just : Plflt, text : PlcharVector) + fun plrandd = c_plrandd : Plflt + fun plreplot = c_plreplot + fun plrgbhls = c_plrgbhls(r : Plflt, g : Plflt, b : Plflt, p_h : PlfltNcScalar, p_l : PlfltNcScalar, p_s : PlfltNcScalar) + fun plschr = c_plschr(def : Plflt, scale : Plflt) + fun plscmap0 = c_plscmap0(r : PlintVector, g : PlintVector, b : PlintVector, ncol0 : Plint) + fun plscmap0a = c_plscmap0a(r : PlintVector, g : PlintVector, b : PlintVector, alpha : PlfltVector, ncol0 : Plint) + fun plscmap0n = c_plscmap0n(ncol0 : Plint) + fun plscmap1 = c_plscmap1(r : PlintVector, g : PlintVector, b : PlintVector, ncol1 : Plint) + fun plscmap1a = c_plscmap1a(r : PlintVector, g : PlintVector, b : PlintVector, alpha : PlfltVector, ncol1 : Plint) + fun plscmap1l = c_plscmap1l(itype : Plbool, npts : Plint, intensity : PlfltVector, coord1 : PlfltVector, coord2 : PlfltVector, coord3 : PlfltVector, alt_hue_path : PlboolVector) + fun plscmap1la = c_plscmap1la(itype : Plbool, npts : Plint, intensity : PlfltVector, coord1 : PlfltVector, coord2 : PlfltVector, coord3 : PlfltVector, alpha : PlfltVector, alt_hue_path : PlboolVector) + fun plscmap1n = c_plscmap1n(ncol1 : Plint) + fun plscmap1_range = c_plscmap1_range(min_color : Plflt, max_color : Plflt) + fun plscol0 = c_plscol0(icol0 : Plint, r : Plint, g : Plint, b : Plint) + fun plscol0a = c_plscol0a(icol0 : Plint, r : Plint, g : Plint, b : Plint, alpha : Plflt) + fun plscolbg = c_plscolbg(r : Plint, g : Plint, b : Plint) + fun plscolbga = c_plscolbga(r : Plint, g : Plint, b : Plint, alpha : Plflt) + fun plscolor = c_plscolor(color : Plint) + fun plscompression = c_plscompression(compression : Plint) + fun plsdev = c_plsdev(devname : PlcharVector) + fun plsdidev = c_plsdidev(mar : Plflt, aspect : Plflt, jx : Plflt, jy : Plflt) + fun plsdimap = c_plsdimap(dimxmin : Plint, dimxmax : Plint, dimymin : Plint, dimymax : Plint, dimxpmm : Plflt, dimypmm : Plflt) + fun plsdiori = c_plsdiori(rot : Plflt) + fun plsdiplt = c_plsdiplt(xmin : Plflt, ymin : Plflt, xmax : Plflt, ymax : Plflt) + fun plsdiplz = c_plsdiplz(xmin : Plflt, ymin : Plflt, xmax : Plflt, ymax : Plflt) + fun plsdrawmode = c_plsdrawmode(mode : Plint) + fun plseed = c_plseed(seed : LibC::UInt) + fun plsesc = c_plsesc(esc : LibC::Char) + fun plsfam = c_plsfam(fam : Plint, num : Plint, bmax : Plint) + fun plsfci = c_plsfci(fci : Plunicode) + fun plsfnam = c_plsfnam(fnam : PlcharVector) + fun plsfont = c_plsfont(family : Plint, style : Plint, weight : Plint) + fun plshade = c_plshade(a : PlfltMatrix, nx : Plint, ny : Plint, defined : PldefinedCallback, xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, shade_min : Plflt, shade_max : Plflt, sh_cmap : Plint, sh_color : Plflt, sh_width : Plflt, min_color : Plint, min_width : Plflt, max_color : Plint, max_width : Plflt, fill : PlfillCallback, rectangular : Plbool, pltr : PltransformCallback, pltr_data : PlPointer) + alias PldefinedCallback = (Plflt, Plflt -> Plint) + alias PlfillCallback = (Plint, PlfltVector, PlfltVector -> Void) + fun plshades = c_plshades(a : PlfltMatrix, nx : Plint, ny : Plint, defined : PldefinedCallback, xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, clevel : PlfltVector, nlevel : Plint, fill_width : Plflt, cont_color : Plint, cont_width : Plflt, fill : PlfillCallback, rectangular : Plbool, pltr : PltransformCallback, pltr_data : PlPointer) + fun plslabelfunc = c_plslabelfunc(label_func : PllabelFuncCallback, label_data : PlPointer) + alias PllabelFuncCallback = (Plint, Plflt, PlcharNcVector, Plint, PlPointer -> Void) + fun plsmaj = c_plsmaj(def : Plflt, scale : Plflt) + fun plsmem = c_plsmem(maxx : Plint, maxy : Plint, plotmem : PlPointer) + fun plsmema = c_plsmema(maxx : Plint, maxy : Plint, plotmem : PlPointer) + fun plsmin = c_plsmin(def : Plflt, scale : Plflt) + fun plsori = c_plsori(ori : Plint) + fun plspage = c_plspage(xp : Plflt, yp : Plflt, xleng : Plint, yleng : Plint, xoff : Plint, yoff : Plint) + fun plspal0 = c_plspal0(filename : PlcharVector) + fun plspal1 = c_plspal1(filename : PlcharVector, interpolate : Plbool) + fun plspause = c_plspause(pause : Plbool) + fun plsstrm = c_plsstrm(strm : Plint) + fun plssub = c_plssub(nx : Plint, ny : Plint) + fun plssym = c_plssym(def : Plflt, scale : Plflt) + fun plstar = c_plstar(nx : Plint, ny : Plint) + fun plstart = c_plstart(devname : PlcharVector, nx : Plint, ny : Plint) + fun plstransform = c_plstransform(coordinate_transform : PltransformCallback, coordinate_transform_data : PlPointer) + fun plstring = c_plstring(n : Plint, x : PlfltVector, y : PlfltVector, string : PlcharVector) + fun plstring3 = c_plstring3(n : Plint, x : PlfltVector, y : PlfltVector, z : PlfltVector, string : PlcharVector) + fun plstripa = c_plstripa(id : Plint, pen : Plint, x : Plflt, y : Plflt) + fun plstripc = c_plstripc(id : PlintNcScalar, xspec : PlcharVector, yspec : PlcharVector, xmin : Plflt, xmax : Plflt, xjump : Plflt, ymin : Plflt, ymax : Plflt, xlpos : Plflt, ylpos : Plflt, y_ascl : Plbool, acc : Plbool, colbox : Plint, collab : Plint, colline : PlintVector, styline : PlintVector, legline : PlcharMatrix, labx : PlcharVector, laby : PlcharVector, labtop : PlcharVector) + fun plstripd = c_plstripd(id : Plint) + fun plimagefr = c_plimagefr(idata : PlfltMatrix, nx : Plint, ny : Plint, xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, zmin : Plflt, zmax : Plflt, valuemin : Plflt, valuemax : Plflt, pltr : PltransformCallback, pltr_data : PlPointer) + fun plimage = c_plimage(idata : PlfltMatrix, nx : Plint, ny : Plint, xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, zmin : Plflt, zmax : Plflt, dxmin : Plflt, dxmax : Plflt, dymin : Plflt, dymax : Plflt) + fun plstyl = c_plstyl(nms : Plint, mark : PlintVector, space : PlintVector) + fun plsurf3d = c_plsurf3d(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint, clevel : PlfltVector, nlevel : Plint) + fun plsurf3dl = c_plsurf3dl(x : PlfltVector, y : PlfltVector, z : PlfltMatrix, nx : Plint, ny : Plint, opt : Plint, clevel : PlfltVector, nlevel : Plint, indexxmin : Plint, indexxmax : Plint, indexymin : PlintVector, indexymax : PlintVector) + fun plsvect = c_plsvect(arrowx : PlfltVector, arrowy : PlfltVector, npts : Plint, fill : Plbool) + fun plsvpa = c_plsvpa(xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt) + fun plsxax = c_plsxax(digmax : Plint, digits : Plint) + fun plsyax = c_plsyax(digmax : Plint, digits : Plint) + fun plsym = c_plsym(n : Plint, x : PlfltVector, y : PlfltVector, code : Plint) + fun plszax = c_plszax(digmax : Plint, digits : Plint) + fun pltext = c_pltext + fun pltimefmt = c_pltimefmt(fmt : PlcharVector) + fun plvasp = c_plvasp(aspect : Plflt) + fun plvect = c_plvect(u : PlfltMatrix, v : PlfltMatrix, nx : Plint, ny : Plint, scale : Plflt, pltr : PltransformCallback, pltr_data : PlPointer) + fun plvpas = c_plvpas(xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, aspect : Plflt) + fun plvpor = c_plvpor(xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt) + fun plvsta = c_plvsta + fun plw3d = c_plw3d(basex : Plflt, basey : Plflt, height : Plflt, xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt, zmin : Plflt, zmax : Plflt, alt : Plflt, az : Plflt) + fun plwidth = c_plwidth(width : Plflt) + fun plwind = c_plwind(xmin : Plflt, xmax : Plflt, ymin : Plflt, ymax : Plflt) + fun plxormod = c_plxormod(mode : Plbool, status : PlboolNcScalar) + alias PlboolNcScalar = Plbool* + fun plsetopt = c_plsetopt(opt : PlcharVector, optarg : PlcharVector) : Plint + fun plparseopts = c_plparseopts(p_argc : LibC::Int*, argv : PlcharNcMatrix, mode : Plint) : Plint + alias PlcharNcMatrix = LibC::Char** +end diff --git a/src/nn/datasets/iris.cr b/src/nn/datasets/iris.cr new file mode 100644 index 00000000..dae37949 --- /dev/null +++ b/src/nn/datasets/iris.cr @@ -0,0 +1,62 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +require "http" + +module Num::NN + extend self + + IRIS_URL = "https://rawgist.apw.app/curran/a08a1080b88344b0c8a7/raw/639388c2cbc2120a14dcf466e85730eb8be498bb/iris.csv" + + def load_iris_dataset + response = HTTP::Client.get(IRIS_URL) + csv = CSV.parse response.body + + features = csv[1...].map &.[...-1] + labels = csv[1...].map &.[-1] + + rng = (0...labels.size).to_a + rng.shuffle! + + features = features.map_with_index do |_, i| + features[rng[i]] + end + + labels = labels.map_with_index do |_, i| + labels[rng[i]] + end + + x_train = features.to_tensor.as_type(Float64) + + label_mapping = { + "setosa" => [0, 0, 1], + "versicolor" => [0, 1, 0], + "virginica" => [1, 0, 0], + } + + mapped = labels.map { |el| label_mapping[el] } + y_train = mapped.to_tensor.as_type(Float64) + + {labels, x_train, y_train} + end +end diff --git a/src/nn/gates/convolution.cr b/src/nn/gates/convolution.cr new file mode 100644 index 00000000..3dbe22cb --- /dev/null +++ b/src/nn/gates/convolution.cr @@ -0,0 +1,62 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::ConvolutionGate(T) < Num::Grad::Gate(T) + getter cached_input : Num::Grad::Variable(T) + getter weight : Num::Grad::Variable(T) + getter bias : Num::Grad::Variable(T) + getter padding : Tuple(Int32, Int32) + getter stride : Tuple(Int32, Int32) + + def initialize( + @cached_input : Num::Grad::Variable(T), + @weight : Num::Grad::Variable(T), + @bias : Num::Grad::Variable(T), + @padding : Tuple(Int32, Int32), + @stride : Tuple(Int32, Int32) + ) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + r0, r1, r2 = Num::NN.conv2d_backward( + @cached_input.value, + @weight.value, + @bias.value, gradient, + @padding, + @stride + ) + + [r0, r1, r2] + end + + def cache(result : Num::Grad::Variable(T), *args) + input, weight, bias, padding, stride = args + + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Conv", self, result, input, weight, bias) + end +end diff --git a/src/nn/gates/elu.cr b/src/nn/gates/elu.cr new file mode 100644 index 00000000..7e6df3c4 --- /dev/null +++ b/src/nn/gates/elu.cr @@ -0,0 +1,41 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::EluGate(T) < Num::Grad::Gate(T) + getter cache : T + + def initialize(@cache : T) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [Num::NN.elu_prime(gradient, @cache)] + end + + def cache(result : Num::Grad::Variable(T), *args : Num::Grad::Variable(T)) + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Elu", self, result, *args) + end +end diff --git a/src/nn/gates/flatten.cr b/src/nn/gates/flatten.cr new file mode 100644 index 00000000..76ad9648 --- /dev/null +++ b/src/nn/gates/flatten.cr @@ -0,0 +1,44 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::FlattenGate(T) < Num::Grad::Gate(T) + getter input : Num::Grad::Variable(T) + getter cached_shape : Array(Int32) + + def initialize(@input : Num::Grad::Variable(T), @cached_shape : Array(Int32)) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [gradient.reshape([gradient.shape[0]] + @cached_shape)] + end + + def cache(result : Num::Grad::Variable(T), *args) + a = args[0] + + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Reshape", self, result, a) + end +end diff --git a/src/nn/gates/leaky_relu.cr b/src/nn/gates/leaky_relu.cr new file mode 100644 index 00000000..1e632a06 --- /dev/null +++ b/src/nn/gates/leaky_relu.cr @@ -0,0 +1,41 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::LeakyReluGate(T) < Num::Grad::Gate(T) + getter cache : T + + def initialize(@cache : T) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [Num::NN.leaky_relu_prime(gradient, @cache)] + end + + def cache(result : Num::Grad::Variable(T), *args : Num::Grad::Variable(T)) + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("LeakyRelu", self, result, *args) + end +end diff --git a/src/frame/frame_slice.cr b/src/nn/gates/linear.cr similarity index 55% rename from src/frame/frame_slice.cr rename to src/nn/gates/linear.cr index 142d4fb0..576c59af 100644 --- a/src/frame/frame_slice.cr +++ b/src/nn/gates/linear.cr @@ -21,44 +21,41 @@ # OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION # WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. -require "./frame" +class Num::NN::LinearGate(T) < Num::Grad::Gate(T) + getter input : Num::Grad::Variable(T) + getter weight : Num::Grad::Variable(T) + getter bias : Num::Grad::Variable(T) -# A `FrameSlice` is a row copy of a `DataFrame`. It owns its -# own data, and cannot be modified. It is the result of a -# row-wise slice of a `DataFrame`. -# -# When used it an operation against a `DataFrame`, it broadcasts -# against the columns of a `DataFrame`, so that reductions will -# always be able to be operated upon -class FrameSlice(T) - getter c : T - - # Initializes a DataFrame from a variadic number of arguments. - # This should only be used by the `DataFrame` class, but remains - # public - def initialize(**args : **T) - @c = args + def initialize(@input : Num::Grad::Variable(T), @weight : Num::Grad::Variable(T), @bias : Num::Grad::Variable(T)) end - # :nodoc: - def to_s(io) - kw = 0 - vw = 0 - @c.each do |k, v| - ks = "#{k}".size - vs = Num::Internal.format(v).size - if ks > kw - kw = ks - end - if vs > vw - vw = vs - end + def backward(payload : Num::Grad::Payload(T)) : Array(T) + grad = payload.variable.grad + + result = [ + grad, + grad, + grad, + ] + + if @input.requires_grad + result[0] = grad.matmul(@weight.value) + end + + if @weight.requires_grad + result[1] = grad.transpose.matmul(@input.value) end - @c.each do |k, v| - io << "#{k}".rjust(kw) - io << " " - io << Num::Internal.format(v).rjust(vw) - io << "\n" + + if @bias.requires_grad + result[2] = grad.sum(axis: 0) end + + result + end + + def cache(result : Num::Grad::Variable(T), *args : Num::Grad::Variable(T)) + input, weight, bias = args + result.grad = T.zeros_like(result.value) + Num::Grad.register("Linear", self, result, input, weight, bias) end end diff --git a/src/nn/gates/mse.cr b/src/nn/gates/mse.cr new file mode 100644 index 00000000..5708dbb1 --- /dev/null +++ b/src/nn/gates/mse.cr @@ -0,0 +1,52 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::MSEGate(T) < Num::Grad::Gate(T) + getter target : T + getter cache : Num::Grad::Variable(T) + + def initialize(@target : T, @cache : Num::Grad::Variable(T)) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + grad = gradient.value + norm = grad * 2 / gradient.size + + output = @cache.value.map(target) do |x, y| + norm * (x - y) + end + + [output] + end + + def cache(result : Num::Grad::Variable(T), *args) + input, target = args + + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("MSE", self, result, input) + end +end diff --git a/src/nn/gates/relu.cr b/src/nn/gates/relu.cr new file mode 100644 index 00000000..8993da17 --- /dev/null +++ b/src/nn/gates/relu.cr @@ -0,0 +1,41 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::ReluGate(T) < Num::Grad::Gate(T) + getter cache : T + + def initialize(@cache : T) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [Num::NN.relu_prime(gradient, @cache)] + end + + def cache(result : Num::Grad::Variable(T), *args : Num::Grad::Variable(T)) + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Relu", self, result, *args) + end +end diff --git a/src/nn/gates/sigmoid.cr b/src/nn/gates/sigmoid.cr new file mode 100644 index 00000000..e107d654 --- /dev/null +++ b/src/nn/gates/sigmoid.cr @@ -0,0 +1,41 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::SigmoidGate(T) < Num::Grad::Gate(T) + getter cache : T + + def initialize(@cache : T) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + [Num::NN.sigmoid_prime(gradient, @cache)] + end + + def cache(result : Num::Grad::Variable(T), *args : Num::Grad::Variable(T)) + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("Sigmoid", self, result, *args) + end +end diff --git a/src/nn/gates/sigmoid_cross_entropy.cr b/src/nn/gates/sigmoid_cross_entropy.cr new file mode 100644 index 00000000..6946daf5 --- /dev/null +++ b/src/nn/gates/sigmoid_cross_entropy.cr @@ -0,0 +1,53 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::SigmoidCrossEntropy(T) < Num::Grad::Gate(T) + getter target : T + getter cache : Num::Grad::Variable(T) + + def initialize(@target : T, @cache : Num::Grad::Variable(T)) + end + + def backward(payload : Num::Grad::Payload(T)) : Array(T) + gradient = payload.variable.grad + + grad = gradient.value + + batch_size = @cache.value.shape[0] + + output = @cache.value.map(@target) do |x, y| + grad * ((1 / (1 + Math.exp(-x))) - y) / batch_size + end + + [output] + end + + def cache(result : Num::Grad::Variable, *args) + a, target = args + + result.grad = T.zeros_like(result.value) + result.requires_grad = true + + Num::Grad.register("SigmoidCrossEntropy", self, result, a) + end +end diff --git a/src/nn/layers/convolution.cr b/src/nn/layers/convolution.cr new file mode 100644 index 00000000..30aa2107 --- /dev/null +++ b/src/nn/layers/convolution.cr @@ -0,0 +1,60 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::ConvolutionalLayer(T) < Num::NN::Layer(T) + getter weights : Num::Grad::Variable(T) + getter bias : Num::Grad::Variable(T) + getter padding : Tuple(Int32, Int32) + getter stride : Tuple(Int32, Int32) = {1, 1} + + def initialize( + context : Num::Grad::Context(T), + in_shape : Array(Int), + num_filters : Int, + kernel_height : Int, + kernel_width : Int, + @padding = {0, 0}, + @stride = {1, 1} + ) + c_in, h_in, w_in = in_shape + w = T.normal([num_filters, c_in, kernel_height, kernel_width]) + b = T.zeros([num_filters, 1, 1]) + @weights = context.variable(w) + @bias = context.variable(w) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = Num::NN.conv2d(input.value, @weights.value, @bias.value, padding, stride) + result = input.context.variable(output) + + if input.is_grad_needed || @weights.is_grad_needed || @bias.is_grad_needed + gate = Num::NN::ConvolutionGate.new(input, @weights, @bias, @padding, @stride) + gate.cache(result, input, @weights, @bias, @padding, @stride) + end + result + end + + def variables : Array(Num::Grad::Variable(T)) + [weights, bias] + end +end diff --git a/src/nn/layers/elu.cr b/src/nn/layers/elu.cr new file mode 100644 index 00000000..ad0829d0 --- /dev/null +++ b/src/nn/layers/elu.cr @@ -0,0 +1,51 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::EluLayer(T) < Num::NN::Layer(T) + def initialize(context : Num::Grad::Context(T)) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = Num::NN.elu(input.value) + result = input.context.variable(output) + + if input.is_grad_needed + gate = Num::NN::EluGate.new(input.value) + gate.cache(result, input) + end + result + end +end + +class Num::Grad::Variable(T) + def elu + output = Num::NN.elu(@value) + result = @context.variable(output) + + if self.is_grad_needed + gate = Num::NN::EluGate.new(@value) + gate.cache(result, self) + end + result + end +end diff --git a/src/nn/layers/flatten.cr b/src/nn/layers/flatten.cr new file mode 100644 index 00000000..ef428da8 --- /dev/null +++ b/src/nn/layers/flatten.cr @@ -0,0 +1,40 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::FlattenLayer(T) < Num::NN::Layer(T) + getter shape : Array(Int32) + + def initialize(context : Num::Grad::Context(T), @shape : Array(Int32)) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = input.value.reshape([input.value.shape[0], -1]) + result = input.context.variable(output) + + if input.is_grad_needed + gate = Num::NN::FlattenGate.new(result, @shape) + gate.cache(result, input) + end + result + end +end diff --git a/src/nn/layers/leaky_relu.cr b/src/nn/layers/leaky_relu.cr new file mode 100644 index 00000000..e07fe356 --- /dev/null +++ b/src/nn/layers/leaky_relu.cr @@ -0,0 +1,51 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::LeakyReluLayer(T) < Num::NN::Layer(T) + def initialize(context : Num::Grad::Context(T)) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = Num::NN.leaky_relu(input.value) + result = input.context.variable(output) + + if input.is_grad_needed + gate = Num::NN::LeakyReluGate.new(input.value) + gate.cache(result, input) + end + result + end +end + +class Num::Grad::Variable(T) + def leaky_relu + output = Num::NN.leaky_relu(@value) + result = @context.variable(output) + + if self.is_grad_needed + gate = Num::NN::LeakyReluGate.new(@value) + gate.cache(result, self) + end + result + end +end diff --git a/src/nn/layers/linear.cr b/src/nn/layers/linear.cr new file mode 100644 index 00000000..9a158bbd --- /dev/null +++ b/src/nn/layers/linear.cr @@ -0,0 +1,49 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::LinearLayer(T) < Num::NN::Layer(T) + getter weights : Num::Grad::Variable(T) + getter bias : Num::Grad::Variable(T) + + def initialize(context : Num::Grad::Context(T), inp_dim : Int, outp_dim : Int) + w = T.rand([outp_dim, inp_dim]) + b = T.zeros([1, outp_dim]) + @weights = context.variable(w) + @bias = context.variable(b) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = input.value.matmul(@weights.value.transpose) + bias.value + result = input.context.variable(output) + + if input.is_grad_needed || @weights.is_grad_needed || @bias.is_grad_needed + gate = Num::NN::LinearGate.new(input, @weights, @bias) + gate.cache(result, input, @weights, @bias) + end + result + end + + def variables : Array(Num::Grad::Variable(T)) + [weights, bias] + end +end diff --git a/src/nn/layers/relu.cr b/src/nn/layers/relu.cr new file mode 100644 index 00000000..1bc33dd4 --- /dev/null +++ b/src/nn/layers/relu.cr @@ -0,0 +1,51 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::ReluLayer(T) < Num::NN::Layer(T) + def initialize(context : Num::Grad::Context(T)) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = Num::NN.relu(input.value) + result = input.context.variable(output) + + if input.is_grad_needed + gate = Num::NN::ReluGate.new(input.value) + gate.cache(result, input) + end + result + end +end + +class Num::Grad::Variable(T) + def relu + output = Num::NN.relu(@value) + result = @context.variable(output) + + if self.is_grad_needed + gate = Num::NN::ReluGate.new(@value) + gate.cache(result, self) + end + result + end +end diff --git a/src/nn/layers/sigmoid.cr b/src/nn/layers/sigmoid.cr new file mode 100644 index 00000000..f46ba9f2 --- /dev/null +++ b/src/nn/layers/sigmoid.cr @@ -0,0 +1,38 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::SigmoidLayer(T) < Num::NN::Layer(T) + def initialize(context : Num::Grad::Context(T)) + end + + def forward(input : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + output = Num::NN.sigmoid(input.value) + result = input.context.variable(output) + + if input.is_grad_needed + gate = Num::NN::SigmoidGate.new(input.value) + gate.cache(result, input) + end + result + end +end diff --git a/src/nn/loss.cr b/src/nn/loss.cr new file mode 100644 index 00000000..66852434 --- /dev/null +++ b/src/nn/loss.cr @@ -0,0 +1,50 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::SigmoidCrossEntropyLoss(T) < Num::NN::Loss(T) + def loss(input : Num::Grad::Variable(T), target : T) : Num::Grad::Variable(T) + output = Num::NN.sigmoid_cross_entropy(input.value, target) + + result = input.context.variable([output]) + + if input.is_grad_needed + gate = Num::NN::SigmoidCrossEntropy(T).new(target, input) + gate.cache(result, input, target) + end + result + end +end + +class Num::NN::MSELoss(T) < Num::NN::Loss(T) + def loss(input : Num::Grad::Variable(T), target : T) : Num::Grad::Variable(T) + output = Num::NN.mse(input.value, target) + + result = input.context.variable(output) + + if input.is_grad_needed + gate = Num::NN::MSEGate(T).new(target, input) + gate.cache(result, input, target) + end + result + end +end diff --git a/src/nn/network.cr b/src/nn/network.cr new file mode 100644 index 00000000..1ac2f515 --- /dev/null +++ b/src/nn/network.cr @@ -0,0 +1,346 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::NetworkInfo(T) + getter layers : Array(Num::NN::Layer(T)) + getter optimizer : Num::NN::Optimizer(T) + getter context : Num::Grad::Context(T) + getter loss : Num::NN::Loss(T) + + # This should always be initialized with an empty + # array of layers that can be tapped and yielded + # by Network creation + def initialize(@context : Num::Grad::Context(T)) + @layers = [] of Num::NN::Layer(T) + @optimizer = Num::NN::Optimizer(T).new + @loss = Num::NN::Loss(T).new + end + + # Add a linear layer to the Network. Since activation functions + # are just treated as additional layers, this simply requires + # the dimensions of the transformation. + # + # Dimensions should be `NUM_FEATURES` x `NUM_OUTPUTS`, so + # if the data set is 100x10, with 200 neurons in the hidden layers, + # the dimensions of the layer would be 10, 100, the 200 will be handled + # by dynamically. + # + # Arguments + # --------- + # *i* : Int + # The number of features in the linear layer + # *j* : Int + # The number of outputs in the linear layer + # + # Examples + # -------- + # ``` + # net = Num::NN::Network.new(ctx) do + # linear 2, 3 + # end + # ``` + def linear(i : Int, j : Int) + @layers << Num::NN::LinearLayer(T).new(@context, i, j) + end + + # Convolution layer for a neural network + # + # Arguments + # --------- + # shape : Array(Int32) + # Shape of the input to the layer + # n : Int32 + # Number of filters to apply + # kh : Int32 + # Filter height + # kw : Int32 + # Filter width + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def conv2d(shape : Array(Int32), n : Int32, kh : Int32, kw : Int32) + @layers << Num::NN::ConvolutionalLayer(T).new(@context, shape, n, kh, kw) + end + + # Adds a Flattening layer to a neural network + # + # Arguments + # --------- + # shape : Array(Int32) + # Shape of input to the layer + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def flatten(shape : Array(Int32)) + @layers << Num::NN::FlattenLayer(T).new(@context, shape) + end + + # Add a ReLU layer to the Network. Activation functions are handled + # the same way as other layers, but do not change the dimensions + # of the input + # + # Arguments + # --------- + # + # Examples + # -------- + # ``` + # net = Num::NN::Network.new(ctx) do + # linear 2, 3 + # relu + # end + # ``` + def relu + @layers << Num::NN::ReluLayer(T).new(@context) + end + + # Adds a Leaky ReLU layer to a network. + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def leaky_relu + @layers << Num::NN::LeakyReluLayer(T).new(@context) + end + + # Adds an ELU layer to the network + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def elu + @layers << Num::NN::EluLayer(T).new(@context) + end + + # Add a Sigmoid layer to the Network. Activation functions are handled + # the same way as other layers, but do not change the dimensions + # of the input + # + # Arguments + # --------- + # + # Examples + # -------- + # ``` + # net = Num::NN::Network.new(ctx) do + # linear 2, 3 + # sigmoid + # end + # ``` + def sigmoid + @layers << Num::NN::SigmoidLayer(T).new(@context) + end + + # Add an SGD optimizer to the Network. + # + # Arguments + # --------- + # *learning_rate* : Float64 + # Learning rate for all layers in the Network + # + # Examples + # -------- + # ``` + # net = Num::NN::Network.new(ctx) do + # linear 2, 3 + # sigmoid + # linear 3, 1 + # sgd 0.7 + # end + # ``` + def sgd(learning_rate : Float64 = 0.01) + @optimizer = Num::NN::SGDOptimizer(T).new(learning_rate) + end + + def adam(*args) + @optimizer = Num::NN::AdamOptimizer(T).new(*args) + end + + # Uses Sigmoid Cross Entropy to compute the loss for + # the Network + # + # Arguments + # --------- + # + # Examples + # -------- + # ``` + # net = Num::NN::Network.new(ctx) do + # linear 2, 3 + # sigmoid + # linear 3, 1 + # sgd 0.7 + # sigmoid_cross_entropy_loss + # end + # ``` + def sigmoid_cross_entropy_loss + @loss = Num::NN::SigmoidCrossEntropyLoss(T).new + end + + # Uses Mean Squared Error to compute the loss for + # the Network + # + # Arguments + # --------- + # + # Examples + # -------- + # ``` + # net = Num::NN::Network.new(ctx) do + # linear 2, 3 + # sigmoid + # linear 3, 1 + # sgd 0.7 + # mse_loss + # end + # ``` + def mse_loss + @loss = Num::NN::MSELoss(T).new + end + + forward_missing_to layers +end + +# A neural network can be defined as a biologically inspired +# computational model that consists of a network architecture +# composed of artificial neurons. This structure contains a set of +# parameters, which can be adjusted to perform specific tasks. +# +# This class is a loose wrapper that primarily provides syntactic +# sugar around moving data through a network -> forward, and propogating +# loss <- backwards +class Num::NN::Network(T) + @layers : Num::NN::NetworkInfo(T) + + # Convenience method to allow for creation of a Network + # with as little code as possible. Taps an instance of + # a LayerArray in order to allow layers to be added to the + # network in a block + # + # Arguments + # --------- + # + # Examples + # -------- + # ``` + # Network(Float32).new do + # layer(2, 3, :tanh) + # layer(3, 1, :sigmoid) + # end + # ``` + def self.new(context : Num::Grad::Context(T), **options) + layers = Num::NN::NetworkInfo(T).new(context) + layers.tap do |instance| + with instance yield + end + layers.optimizer.build_params(layers.layers) + new(layers, **options) + end + + # :nodoc: + private def initialize(@layers : Num::NN::NetworkInfo(T), **options) + end + + # Propogates an input through a network, returning + # the final prediction from the network + # + # Arguments + # --------- + # *train* : Tensor(T) + # Training input data + # + # Examples + # -------- + # ``` + # a = [[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]].to_tensor + # net.forward(a) + # ``` + def forward(train : Num::Grad::Variable(T)) : Num::Grad::Variable(T) + @layers.each do |layer| + train = layer.forward(train) + end + train + end + + # Uses the Network's loss function to calculate the loss + # based on the final output from the Network, as well + # as the target output + # + # Arguments + # --------- + # *output* : Num::Grad::Variable(T) + # Prediction by the network + # *target* : T + # Tensor containing ground truth values + # + # Examples + # -------- + # ``` + # epochs.times do |epoch| + # y_pred = net.forward(x) + # loss = net.loss(y_pred, y_actual) + # end + # ``` + def loss(output : Num::Grad::Variable(T), target : T) + @layers.loss.loss(output, target) + end + + # Return the Network's optimizer to allow updating + # the weights and biases of the network + # + # Arguments + # --------- + # + # Examples + # -------- + # ``` + # epochs.times do |epoch| + # y_pred = net.forward(x) + # loss = net.loss(y_pred, y_actual) + # net.optimizer.update + # end + # ``` + def optimizer + @layers.optimizer + end +end diff --git a/src/nn/optimizer.cr b/src/nn/optimizer.cr new file mode 100644 index 00000000..7bb8a4d5 --- /dev/null +++ b/src/nn/optimizer.cr @@ -0,0 +1,168 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::SGDOptimizer(T) < Num::NN::Optimizer(T) + getter params : Array(Num::Grad::Variable(T)) + getter learning_rate : Float64 + + def initialize(@learning_rate : Float64) + @params = [] of Num::Grad::Variable(T) + end + + def build_params(l : Array(Layer(T))) + l.each do |layer| + layer.variables.each do |v| + @params << v + end + end + end + + def update + @params.each do |v| + if v.requires_grad + v.value.map!(v.grad) do |x, y| + x - @learning_rate * y + end + v.grad = T.zeros_like(v.value) + end + end + end +end + +class Num::NN::SGDMomentumOptimizer(T) < Num::NN::Optimizer(T) + getter params : Array(Num::Grad::Variable(T)) + getter learning_rate : Float64 + getter momentum : Float64 + getter moments : Array(T) + getter decay : Float64 + getter nesterov : Bool + + def initialize( + @learning_rate : Float64, + @momentum : Float64 = 0.0, + @decay : Float64 = 0.0, + @nesterov : Bool = false + ) + @params = [] of Num::Grad::Variable(T) + @moments = [] of T + end + + def build_params(l : Array(Layer(T))) + l.each do |layer| + layer.variables.each do |v| + @params << v + @moments << T.zeros_like(v.grad) + end + end + end + + def update + @learning_rate *= 1 / (@decay + 1) + @decay += @decay + + @params.size.times do |i| + v = @params[i] + + if v.requires_grad + @moments[i].map!(v.grad) do |x, y| + @momentum * x - @learning_rate * y + end + + if @nesterov + v.value.map!(v.grad, @moments[i]) do |x, y, z| + x - @learning_rate * y + @momentum * z + end + else + v.value.map!(@moments[i]) do |x, y| + x + y + end + end + + v.grad = T.zeros_like(v.value) + end + end + end +end + +class Num::NN::AdamOptimizer(T) < Num::NN::Optimizer(T) + getter params : Array(Num::Grad::Variable(T)) + getter learning_rate : Float64 + getter beta1 : Float64 + getter beta2 : Float64 + getter beta1_t : Float64 + getter beta2_t : Float64 + getter first_moments : Array(T) + getter second_moments : Array(T) + getter epsilon : Float64 + + def initialize( + @learning_rate : Float64 = 0.001, + @beta1 : Float64 = 0.9, + @beta2 : Float64 = 0.999, + @epsilon : Float64 = 1e-8 + ) + @params = [] of Num::Grad::Variable(T) + @beta1_t = @beta1 + @beta2_t = @beta2 + + @first_moments = [] of T + @second_moments = [] of T + end + + def build_params(l : Array(Layer(T))) + l.each do |layer| + layer.variables.each do |v| + @params << v + @first_moments << T.zeros_like(v.grad) + @second_moments << T.zeros_like(v.grad) + end + end + end + + def update + lr_t = @learning_rate * Math.sqrt(1 - @beta2_t) / (1 - @beta1_t) + + @beta1_t *= @beta1 + @beta2_t *= @beta2 + + @params.size.times do |i| + v = @params[i] + + if v.requires_grad + @first_moments[i].map!(v.grad) do |x, y| + @beta1 * x + (1 - @beta1) * y + end + + @second_moments[i].map!(v.grad) do |x, y| + @beta2 * x + (1 - @beta2) * y * y + end + + v.value.map!(@first_moments[i], @second_moments[i]) do |x, y, z| + x - lr_t * y / (Math.sqrt(z) + @epsilon) + end + + v.grad = T.zeros_like(v.value) + end + end + end +end diff --git a/src/nn/primitives/activation.cr b/src/nn/primitives/activation.cr new file mode 100644 index 00000000..96fac78c --- /dev/null +++ b/src/nn/primitives/activation.cr @@ -0,0 +1,316 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +module Num::NN + extend self + + # Tanh squashes a real-valued number to the range [-1, 1]. It’s non-linear. + # But unlike Sigmoid, its output is zero-centered. Therefore, in practice + # the tanh non-linearity is always preferred to the sigmoid nonlinearity. + # + # Arguments + # --------- + # *x* : Tensor + # Tensor to activate + # + # Examples + # -------- + # ``` + # a = [0.1, 0.34, 0.65].to_tensor + # Num::NN.tanh(a) # => [0.099668, 0.327477, 0.57167 ] + # ``` + def tanh(x : Tensor(U)) : Tensor(U) forall U + Num.tanh(x) + end + + # :ditto: + def tanh!(x : Tensor(U)) forall U + Num.tanh!(x) + end + + # Derivative of the Tanh function + # + # Arguments + # --------- + # *x* : Tensor + # Tensor to derive + # + # Examples + # -------- + # ``` + # a = [0.1, 0.34, 0.65].to_tensor + # Num::NN.d_tanh(a) # => [0.990066, 0.892759, 0.673193] + # ``` + def tanh_prime(gradient : Tensor(U), cached : Tensor(U)) forall U + cached.map(gradient) do |x, y| + y * (1 - x * x) + end + end + + # Sigmoid takes a real value as input and outputs another value + # between 0 and 1. It’s easy to work with and has all the + # nice properties of activation functions: it’s non-linear, + # continuously differentiable, monotonic, and has a fixed + # output range. + # + # Arguments + # --------- + # *x* : Tensor + # Tensor to activate + # + # Examples + # -------- + # ``` + # a = [0.1, 0.34, 0.65].to_tensor + # puts Num::NN.sigmoid(a) # => [0.524979, 0.584191, 0.65701 ] + # ``` + def sigmoid(x) + x.map do |i| + 1 / (1 + Math.exp(-i)) + end + end + + # :ditto: + def sigmoid!(x : Tensor(U)) : Tensor(U) forall U + x.map! do |i| + 1 / (1 + Math.exp(-i)) + end + end + + # Derivative of the Sigmoid function + # + # Arguments + # --------- + # *x* : Tensor + # Tensor to derive + # + # Examples + # -------- + # ``` + # a = [0.1, 0.34, 0.65].to_tensor + # puts Num::NN.d_sigmoid(a) # => [0.249376, 0.242912, 0.225348] + # ``` + def sigmoid_prime(gradient : Tensor(U), cached : Tensor(U)) : Tensor(U) forall U + cached.map(gradient) do |x, y| + x * (1 - x) * y + end + end + + # ReLU activation function + # + # Arguments + # --------- + # x : Tensor(U) + # Argument to activate + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def relu(x : Tensor(U)) : Tensor(U) forall U + Num.max(U.new(0), x) + end + + # :ditto: + def relu!(x : Tensor(U)) : Tensor(U) forall U + Num.max!(U.new(0), x) + end + + # Derivative of the ReLU activation function + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def relu_prime(gradient : Tensor(U), cached : Tensor(U)) : Tensor(U) forall U + cached.map(gradient) do |x, y| + if x <= 0 + U.new(0) + else + y + end + end + end + + # Leaky ReLU activation function + # + # Arguments + # --------- + # x : Tensor(U) + # Argument to activate + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def leaky_relu(x : Tensor(U)) : Tensor(U) forall U + x.map do |i| + i > 0 ? i : U.new(i * 0.01) + end + end + + # :ditto: + def leaky_relu!(x : Tensor(U)) forall U + x.map! do |i| + i > 0 ? i : i * 0.01 + end + end + + # Brief description of leakyreluprime + # + # Arguments + # --------- + # gradient : Tensor(U) + # Gradient value + # cached : Tensor(U) + # Stored value + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def leaky_relu_prime(gradient : Tensor(U), cached : Tensor(U)) : Tensor(U) forall U + cached.map(gradient) do |x, y| + if x <= 0 + y * 0.01 + else + y + end + end + end + + # Exponential linear unit activation + # + # Arguments + # --------- + # x : Tensor(U) + # Argument to activate + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def elu(x : Tensor(U), alpha = 0.01) : Tensor(U) forall U + x.map do |i| + if i > 0 + i + else + U.new(alpha * (Math::E ** i - 1)) + end + end + end + + # :ditto: + def elu!(x : Tensor(U), alpha = 0.01) : Tensor(U) forall U + x.map! do |i| + if i > 0 + i + else + alpha * (Math::E ** i - 1) + end + end + end + + # Derivative of the ELU activation + # + # Arguments + # --------- + # gradient : Tensor(U) + # Gradient value + # cached : Tensor(U) + # Stored value + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def elu_prime(gradient : Tensor(U), cached : Tensor(U)) : Tensor(U) forall U + cached.map(gradient) do |x, y| + if x <= 0 + Math.exp(y) + else + U.new(1) + end + end + end + + # Sigmoid cross entropy loss + # + # Arguments + # --------- + # input : Tensor(U) + # Predicted values + # target : Tensor(U) + # Truth values + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def sigmoid_cross_entropy(input : Tensor(U), target : Tensor(U)) : U forall U + batch_size = input.shape[0] + result = input.map(target) do |x, y| + -y * x + Math.max(x, U.new(0)) + Math.log1p(Math.exp(-x.abs)) + end + result.sum / U.new(batch_size) + end + + # Mean squared error loss + # + # Arguments + # --------- + # input : Tensor(U) + # Predicted values + # target : Tensor(U) + # Truth values + # + # Returns + # ------- + # Tensor(U) + # + # Examples + # -------- + def mse(input : Tensor(U), target : Tensor(U)) : Tensor(U) forall U + result = input.map(target) do |i, j| + (i - j) ** 2 + end + [U.new(result.mean)].to_tensor + end +end diff --git a/src/nn/primitives/layer.cr b/src/nn/primitives/layer.cr new file mode 100644 index 00000000..88477f3b --- /dev/null +++ b/src/nn/primitives/layer.cr @@ -0,0 +1,30 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +abstract class Num::NN::Layer(T) + abstract def forward(input : Num::Grad::Variable(T)) + + def variables : Array(Num::Grad::Variable(T)) + [] of Num::Grad::Variable(T) + end +end diff --git a/src/nn/primitives/loss.cr b/src/nn/primitives/loss.cr new file mode 100644 index 00000000..2fc855f4 --- /dev/null +++ b/src/nn/primitives/loss.cr @@ -0,0 +1,28 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::Loss(T) + def loss(input : Num::Grad::Variable(T), target : T) : Num::Grad::Variable(T) + input + end +end diff --git a/src/nn/primitives/nnpack_conv.cr b/src/nn/primitives/nnpack_conv.cr new file mode 100644 index 00000000..f5d684da --- /dev/null +++ b/src/nn/primitives/nnpack_conv.cr @@ -0,0 +1,144 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +module Num::NN + def conv2d( + input : Tensor(Float32), + weight : Tensor(Float32), + bias : Tensor(Float32), + padding : Tuple(Int, Int), + stride : Tuple(Int, Int) = {1, 1} + ) + batch_size = input.shape[0] + output_channels = weight.shape[-4] + output_height = (2 * padding[0] + input.shape[-2]) - (weight.shape[-2] - 1) + output_width = (2 * padding[1] + input.shape[-1]) - (weight.shape[-1] - 1) + + input = input.dup(Num::RowMajor) + weight = weight.dup(Num::RowMajor) + + result = Tensor(Float32).new([input.shape[0], output_channels, output_height, output_width]) + + status = LibNNPACK.nnp_convolution_output( + LibNNPACK::NNPConvolutionAlgorithm::AUTO, + batch_size.to_u64, + input.shape[-3].to_u64, + output_channels.to_u64, + LibNNPACK::NNPSize.new(height: input.shape[-2], width: input.shape[-1]), + LibNNPACK::NNPPadding.new(top: padding[0], bottom: padding[0], left: padding[1], right: padding[1]), + LibNNPACK::NNPSize.new(height: weight.shape[-2], width: weight.shape[-1]), + input.to_unsafe, + weight.to_unsafe, + bias.to_unsafe, + result.to_unsafe, + nil, + nil, + LibNNPACK::NNPActivation::IDENTITY, + nil, + nil, + out profile, + ) + + unless status == LibNNPACK::NNPStatus::SUCCESS + raise Exception.new "NNPACK failed with #{status}. Did you run with the -Dnnpack flag?" + end + + result + end + + def conv2d_backward( + input : Tensor(Float32), + weight : Tensor(Float32), + bias : Tensor(Float32), + grad_output : Tensor(Float32), + padding : Tuple(Int, Int), + stride : Tuple(Int, Int) = {1, 1} + ) + batch_size = input.shape[0] + input_channels = input.shape[-3] + output_channels = weight.shape[-4] + output_height = (2 * padding[0] + input.shape[-2]) - (weight.shape[-2] - 1) + output_width = (2 * padding[1] + input.shape[-1]) - (weight.shape[-1] - 1) + + nninput_size = LibNNPACK::NNPSize.new(height: input.shape[-2], width: input.shape[-1]) + nnpadding = LibNNPACK::NNPPadding.new(top: padding[0], bottom: padding[0], left: padding[1], right: padding[1]) + nnkernel_size = LibNNPACK::NNPSize.new(height: weight.shape[-2], width: weight.shape[-1]) + + grad_input = Tensor(Float32).zeros(input.shape) + + status = LibNNPACK.nnp_convolution_input_gradient( + LibNNPACK::NNPConvolutionAlgorithm::AUTO, + batch_size, + input_channels, + output_channels, + nninput_size, + nnpadding, + nnkernel_size, + grad_output.to_unsafe, + weight.to_unsafe, + input.to_unsafe, + nil, + nil, + LibNNPACK::NNPActivation::IDENTITY, + nil, + nil, + out input_profile + ) + + unless status == LibNNPACK::NNPStatus::SUCCESS + raise Exception.new "NNPACK failed with #{status}. Did you run with the -Dnnpack flag?" + end + + grad_weight = Tensor(Float32).zeros(input.shape) + + status = LibNNPACK.nnp_convolution_kernel_gradient( + LibNNPACK::NNPConvolutionAlgorithm::AUTO, + batch_size, + input_channels, + output_channels, + nninput_size, + nnpadding, + nnkernel_size, + input.to_unsafe, + grad_output.to_unsafe, + grad_weight.to_unsafe, + nil, + nil, + LibNNPACK::NNPActivation::IDENTITY, + nil, + nil, + out weight_profile + ) + + unless status == LibNNPACK::NNPStatus::SUCCESS + raise Exception.new "NNPACK failed with #{status}. Did you run with the -Dnnpack flag?" + end + + grad_bias = bias + if bias.rank == 3 + grad_bias = grad_bias.sum(3).sum(2).sum(0).reshape(bias.shape) + end + + {grad_input, grad_weight, grad_bias} + end +end diff --git a/src/nn/primitives/optimizer.cr b/src/nn/primitives/optimizer.cr new file mode 100644 index 00000000..1bad67aa --- /dev/null +++ b/src/nn/primitives/optimizer.cr @@ -0,0 +1,37 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::NN::Optimizer(T) + getter params : Array(Num::Grad::Variable(T)) + getter learning_rate : Float64 = 0.01 + + def initialize(@learning_rate : Float64 = 0.01) + @params = [] of Num::Grad::Variable(T) + end + + def build_params(l : Array(Layer(T))) + end + + def update + end +end diff --git a/src/plot/figures/figure.cr b/src/plot/figures/figure.cr new file mode 100644 index 00000000..e800e507 --- /dev/null +++ b/src/plot/figures/figure.cr @@ -0,0 +1,32 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +# The base class that all other figures should inherit from. +# +# A Figure needs to be able to handle it's own plotting, +# as well as be able to return it's bounds to the overall +# plot +abstract class Num::Plot::Figure + abstract def plot + abstract def update_bounds(bounds : Num::Plot::Bounds) : Num::Plot::Bounds +end diff --git a/src/plot/figures/line.cr b/src/plot/figures/line.cr new file mode 100644 index 00000000..7e49244f --- /dev/null +++ b/src/plot/figures/line.cr @@ -0,0 +1,42 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +require "./xy" + +class Num::Plot::LinePlot < Num::Plot::XYPlot + # Plots a line plot on a figure + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def plot + super + LibPlplot.plline(@size, @x.to_unsafe, @y.to_unsafe) + end +end diff --git a/src/plot/figures/scatter.cr b/src/plot/figures/scatter.cr new file mode 100644 index 00000000..a42c6a63 --- /dev/null +++ b/src/plot/figures/scatter.cr @@ -0,0 +1,61 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::Plot::Scatter < Num::Plot::XYPlot + @code : Int32 = 0 + + # Initializes a Scatter plot + # + # Arguments + # --------- + # x + # Tensor like x-axis argument to plot + # y + # Tensor like y-axis argument to plot + # @color : Int32? = nil + # Color code to use + # @code = 0 + # Symbol code to use + # + # Examples + # -------- + def initialize(x, y, @color : Int32? = nil, @code = 0) + super x, y, @color + end + + # Plots a scatter plot + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def plot + super + LibPlplot.plpoin(@size, @x.to_unsafe, @y.to_unsafe, @code) + end +end diff --git a/src/plot/figures/xy.cr b/src/plot/figures/xy.cr new file mode 100644 index 00000000..72653bcf --- /dev/null +++ b/src/plot/figures/xy.cr @@ -0,0 +1,85 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::Plot::XYPlot < Num::Plot::Figure + @x : Tensor(Float64) + @y : Tensor(Float64) + @size : Int32 + @color : Int32? = nil + + # Initializes a basic XY plot + # + # Arguments + # --------- + # x + # Tensor-like x-axis argument + # y + # Tensor-like y-axis argument + # @color : Int32? = nil + # Color code for the plot + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def initialize(x, y, @color : Int32? = nil) + @x = x.to_tensor.as_type(Float64) + @y = y.to_tensor.as_type(Float64) + + unless @x.size == @y.size + raise Num::Internal::ShapeError.new("Inputs must be the same size") + end + + @size = @x.size + end + + # Base plotting method, sets color to a default value + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def plot + unless @color.nil? + LibPlplot.plcol0(@color.unsafe_as(Int32)) + end + end + + def update_bounds(bounds : Num::Plot::Bounds) : Num::Plot::Bounds + x_min, x_max = @x.min, @x.max + y_min, y_max = @y.min, @y.max + + bounds.x_min = {bounds.x_min, x_min}.min + bounds.x_max = {bounds.x_max, x_max}.max + bounds.y_min = {bounds.y_min, y_min}.min + bounds.y_max = {bounds.y_max, y_max}.max + bounds + end +end diff --git a/src/plot/internal/options.cr b/src/plot/internal/options.cr new file mode 100644 index 00000000..f290e287 --- /dev/null +++ b/src/plot/internal/options.cr @@ -0,0 +1,132 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +require "../figures/*" + +class Num::Plot::Options + getter bounds : Num::Plot::Bounds = Num::Plot::Bounds.new + getter figures : Array(Num::Plot::Figure) = [] of Num::Plot::Figure + + private macro tap_prop(name, dtype, default) + @{{name}} : {{dtype}} = {{ default }} + + def {{name}} + @{{name}} + end + + def {{name}}(val : {{dtype}}) + @{{name}} = val + end + end + + tap_prop term, Symbol?, :qtwidget + tap_prop palette, Symbol, :alternate + tap_prop x_label, String, "" + tap_prop y_label, String, "" + tap_prop label, String, "" + + def initialize + end + + # Plots a line plot + # + # Arguments + # --------- + # x + # Tensor-like x-axis data + # y + # Tensor-like y-axis data + # color = nil + # Color code for data + # + # Returns + # ------- + # nil + # + # Examples + # -------- + # ``` + # Num::Plot::Plot.plot do + # line [1, 2, 3], [4, 5, 3] + # end + # ``` + def line(x, y, color = nil) + plt = Num::Plot::LinePlot.new(x, y, color) + @figures << plt + @bounds = plt.update_bounds(@bounds) + end + + # Plots a scatter plot + # + # Arguments + # --------- + # x + # Tensor like x-axis data + # y + # Tensor-like y-axis data + # color = nil + # Color code for plot + # code : Int32 = 1 + # Symbol to use for points + # + # Returns + # ------- + # nil + # + # Examples + # -------- + # ``` + # Num::Plot::Plot.plot do + # scatter [3, 4, 2], [1, 7, 8], 3, 17 + # end + # ``` + def scatter(x, y, color = nil, code : Int32 = 1) + plt = Num::Plot::Scatter.new(x, y, color, code) + @figures << plt + @bounds = plt.update_bounds(@bounds) + end + + # Plots a custom plot. As long as a user defines a class + # that inherits from Figure, this method will ensure it is + # added to the Plot. This allows users to add custom plots + # using the block syntax + # + # Arguments + # --------- + # cls : U.class + # Class to use for plotting + # *args + # Arguments to create a plot + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def custom(cls : U.class, *args) forall U + plt = U.new(*args) + @figures << plt + @bounds = plt.update_bounds(@bounds) + end +end diff --git a/src/plot/internal/palette.cr b/src/plot/internal/palette.cr new file mode 100644 index 00000000..38af2795 --- /dev/null +++ b/src/plot/internal/palette.cr @@ -0,0 +1,63 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +module Num::Plot + # Sets the color index for cmap0 (see the section called “Color Map0”). + # + # 0 black (default background) + # 1 red (default foreground) + # 2 yellow + # 3 green + # 4 aquamarine + # 5 pink + # 6 wheat + # 7 grey + # 8 brown + # 9 blue + # 10 BlueViolet + # 11 cyan + # 12 turquoise + # 13 magenta + # 14 salmon + # 15 white + # Use plscmap0 to change the entire cmap0 color palette and plscol0 to + # change an individual color in the cmap0 color palette. + enum Color + BLACK = 0 + RED = 1 + YELLOW = 2 + GREEN = 3 + AQUAMARINE = 4 + PINK = 5 + WHEAT = 6 + GREY = 7 + BROWN = 8 + BLUE = 9 + BLUE_VIOLET = 10 + CYAN = 11 + TURQUOISE = 12 + MAGENTA = 13 + SALMON = 14 + WHITE = 15 + end +end diff --git a/src/plot/internal/primitives.cr b/src/plot/internal/primitives.cr new file mode 100644 index 00000000..734a8e28 --- /dev/null +++ b/src/plot/internal/primitives.cr @@ -0,0 +1,38 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +struct Num::Plot::Bounds + property x_min : Float64 = 0.0 + property x_max : Float64 = 0.0 + property y_min : Float64 = 0.0 + property y_max : Float64 = 0.0 +end + +module Num::Plot + COLOR_MAPS = { + default: "cmap0_default.pal", + alternate: "cmap0_alternate.pal", + black_white: "cmap0_black_on_white.pal", + white_bg: "cmap0_white_bg.pal", + } +end diff --git a/src/plot/plot.cr b/src/plot/plot.cr new file mode 100644 index 00000000..ca1f0a51 --- /dev/null +++ b/src/plot/plot.cr @@ -0,0 +1,62 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Num::Plot::Plot + # Class method to faciliate plotting of generic plots + # + # Arguments + # --------- + # + # Returns + # ------- + # nil + # + # Examples + # -------- + def self.plot + plotter = Num::Plot::Options.new + plotter.tap do |instance| + with instance yield + end + + if !plotter.term.nil? + LibPlplot.plsdev(plotter.term.to_s) + end + + palette = Num::Plot::COLOR_MAPS[plotter.palette] + + LibPlplot.plspal0(palette) + LibPlplot.plinit + + b = plotter.bounds + LibPlplot.plenv(b.x_min, b.x_max, b.y_min, b.y_max, 0, 0) + LibPlplot.pllab(plotter.x_label, plotter.y_label, plotter.label) + + plotter.figures.each_with_index do |fig, i| + LibPlplot.plcol0(i + 1) + fig.plot + end + + LibPlplot.plend + end +end diff --git a/src/tensor/creation.cr b/src/tensor/creation.cr index 646b23c9..b81e813f 100644 --- a/src/tensor/creation.cr +++ b/src/tensor/creation.cr @@ -25,35 +25,6 @@ require "./tensor" require "./internal/random" class Tensor(T) - # Creates a `Tensor` sampled from a provided range, with a given - # shape. - # - # The generic types of the `Tensor` are inferred from the endpoints - # of the range - # - # Arguments - # --------- - # *r* : Range(U, U) - # Range of values to sample between - # *shape* : Array(Int) - # Shape of returned `Tensor` - # - # Examples - # -------- - # ``` - # Num::Rand.set_seed(0) - # t = Tensor.random(0...10, [2, 2]) - # t - # - # # [[8, 4], - # # [7, 4]] - # ``` - def self.random(r : Range(U, U), shape : Array(Int)) : Tensor(U) forall U - self.new(shape) do - Num::Rand.generator.rand(r) - end - end - # Creates a `Tensor` of a provided shape, filled with 0. The generic type # must be specified. # diff --git a/src/tensor/extensions/number.cr b/src/tensor/extensions/number.cr new file mode 100644 index 00000000..22701fed --- /dev/null +++ b/src/tensor/extensions/number.cr @@ -0,0 +1,49 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +struct Number + macro add_operator(name, operator) + def {{operator.id}}(other : Tensor) + Num.{{name}}(self, other) + end + end + + add_operator add, :+ + add_operator subtract, :- + add_operator multiply, :* + add_operator divide, :/ + add_operator floordiv, :// + add_operator power, :** + add_operator modulo, :% + add_operator left_shift, :<< + add_operator right_shift, :>> + add_operator bitwise_and, :& + add_operator bitwise_or, :| + add_operator bitwise_xor, :^ + add_operator equal, :== + add_operator not_equal, :!= + add_operator greater, :> + add_operator greater_equal, :>= + add_operator less, :< + add_operator less_equal, :<= +end diff --git a/src/tensor/internal/random.cr b/src/tensor/internal/random.cr index 0fb489bc..22d2cfd2 100644 --- a/src/tensor/internal/random.cr +++ b/src/tensor/internal/random.cr @@ -21,10 +21,14 @@ # OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION # WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. +require "alea" + class Num::Rand - class_getter generator = Random.new + class_getter generator = Alea::Random.new + class_getter stdlib_generator = Random.new def self.set_seed(seed) - @@generator = Random.new(seed) + @@generator = Alea::Random.new(seed) + @@stdlib_generator = Random.new(seed) end end diff --git a/src/tensor/internal/yield_iterators.cr b/src/tensor/internal/yield_iterators.cr new file mode 100644 index 00000000..ab3d5ea9 --- /dev/null +++ b/src/tensor/internal/yield_iterators.cr @@ -0,0 +1,137 @@ +require "../tensor" + +macro init_strided_iteration(coord, backstrides, t_shape, t_strides, t_rank) + {{ coord.id }} = Pointer(Int32).malloc({{ t_rank }}, 0) + {{ backstrides.id }} = Pointer(Int32).malloc({{ t_rank }}) + {{ t_rank }}.times do |i| + {{ backstrides.id }}[i] = {{ t_strides }}[i] * ({{ t_shape }}[i] - 1) + end +end + +macro advance_strided_iteration(coord, backstrides, t_shape, t_strides, t_rank, iter_pos) + ({{ t_rank }} - 1).step(to: 0, by: -1) do |k| + if {{ coord.id }}[k] < {{ t_shape }}[k] - 1 + {{ coord.id }}[k] += 1 + {{ iter_pos }} += {{ t_strides }}[k] + break + else + {{ coord.id }}[k] = 0 + {{ iter_pos }} -= {{ backstrides.id }}[k] + end + end +end + +@[AlwaysInline] +def strided_iteration(t : Tensor) + data = t.to_unsafe + if t.flags.contiguous? + t.size.times do |i| + yield i, data + data += 1 + end + else + t_shape, t_strides, t_rank = t.iter_attrs + init_strided_iteration(:coord, :backstrides, t_shape, t_strides, t_rank) + t.size.times do |i| + yield i, data + advance_strided_iteration(:coord, :backstrides, t_shape, t_strides, t_rank, data) + end + end +end + +@[AlwaysInline] +def dual_strided_iteration(t1 : Tensor, t2 : Tensor) + n = t1.size + + t1_contiguous = t1.flags.contiguous? + t2_contiguous = t2.flags.contiguous? + + t1data = t1.to_unsafe + t2data = t2.to_unsafe + + t1_shape, t1_strides, t1_rank = t1.iter_attrs + t2_shape, t2_strides, t2_rank = t2.iter_attrs + + if t1_contiguous && t2_contiguous + n.times do |i| + yield i, t1data, t2data + t1data += 1 + t2data += 1 + end + elsif t1_contiguous + init_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank) + n.times do |i| + yield i, t1data, t2data + t1data += 1 + advance_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank, t2data) + end + elsif t2_contiguous + init_strided_iteration(:t1_coord, :t1_backstrides, t1_shape, t1_strides, t1_rank) + n.times do |i| + yield i, t1data, t2data + advance_strided_iteration(:t1_coord, :t1_backstrides, t1_shape, t1_strides, t1_rank, t1data) + t2data += 1 + end + else + init_strided_iteration(:t1_coord, :t1_backstrides, t1_shape, t1_strides, t1_rank) + init_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank) + n.times do |i| + yield i, t1data, t2data + advance_strided_iteration(:t1_coord, :t1_backstrides, t1_shape, t1_strides, t1_rank, t1data) + advance_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank, t2data) + end + end +end + +@[AlwaysInline] +def tri_strided_iteration(t1 : Tensor, t2 : Tensor, t3 : Tensor) + n = t1.size + + t1_contiguous = t1.flags.contiguous? + t2_contiguous = t2.flags.contiguous? + t3_contiguous = t3.flags.contiguous? + + t1data = t1.to_unsafe + t2data = t2.to_unsafe + t3data = t3.to_unsafe + + t1_shape, t1_strides, t1_rank = t1.iter_attrs + t2_shape, t2_strides, t2_rank = t2.iter_attrs + t3_shape, t3_strides, t3_rank = t3.iter_attrs + + if t1_contiguous && t2_contiguous && t3_contiguous + n.times do |i| + yield i, t1data, t2data, t3data + t1data += 1 + t2data += 1 + t3data += 1 + end + elsif t1_contiguous && t2_contiguous + init_strided_iteration(:t3_coord, :t3_backstrides, t3_shape, t3_strides, t3_rank) + n.times do |i| + yield i, t1data, t2data, t3data + t1data += 1 + t2data += 1 + advance_strided_iteration(:t3_coord, :t3_backstrides, t3_shape, t3_strides, t3_rank, t3data) + end + elsif t1_contiguous + init_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank) + init_strided_iteration(:t3_coord, :t3_backstrides, t3_shape, t3_strides, t3_rank) + n.times do |i| + yield i, t1data, t2data, t3data + t1data += 1 + advance_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank, t2data) + advance_strided_iteration(:t3_coord, :t3_backstrides, t3_shape, t3_strides, t3_rank, t3data) + end + else + init_strided_iteration(:t1_coord, :t1_backstrides, t1_shape, t1_strides, t1_rank) + init_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank) + init_strided_iteration(:t3_coord, :t3_backstrides, t3_shape, t3_strides, t3_rank) + n.times do |i| + yield i, t1data, t2data, t3data + advance_strided_iteration(:t1_coord, :t1_backstrides, t1_shape, t1_strides, t1_rank, t1data) + advance_strided_iteration(:t2_coord, :t2_backstrides, t2_shape, t2_strides, t2_rank, t2data) + advance_strided_iteration(:t3_coord, :t3_backstrides, t3_shape, t3_strides, t3_rank, t3data) + end + end +end diff --git a/src/tensor/linalg.cr b/src/tensor/linalg.cr index 7f1ae04e..2fb080fc 100644 --- a/src/tensor/linalg.cr +++ b/src/tensor/linalg.cr @@ -537,8 +537,12 @@ class Tensor(T) self.is_matrix other.is_matrix + unless self.shape[1] == other.shape[0] + raise Num::Internal::ShapeError.new("Invalid shapes for matrix multiplication: #{@shape}, #{other.shape}") + end + a = @flags.contiguous? || @flags.fortran? ? self : self.dup(Num::RowMajor) - b = other.flags.contiguous? || flags.fortran? ? other : other.dup(Num::RowMajor) + b = other.flags.contiguous? || other.flags.fortran? ? other : other.dup(Num::RowMajor) m = a.shape[0] n = b.shape[1] k = a.shape[1] @@ -663,6 +667,214 @@ class Tensor(T) self.tensordot(b, [[axes_a], [axes_b]]) end + # Compute the matrix exponential using Pade approximation. + # + # Arguments + # --------- + # *self* + # Matrix of which to compute the exponential + # + # Examples + # -------- + # ``` + # a = [[1.0, 2.0], [-1.0, 3.0]].to_tensor * Complex.new(0, 1) + # puts a.expm + # + # # [[0.426459+1.89218j , -2.13721+-0.978113j], + # # [1.06861+0.489056j , -1.71076+0.914063j ]] + # ``` + def expm + self.is_matrix + + a_l1 = self.norm(order: '1') + n_squarings = 0 + + {% if T == Float64 || T == Complex %} + if a_l1 < 1.495585217958292e-002 + u, v = self.pade_three + elsif a_l1 < 2.539398330063230e-001 + u, v = self.pade_five + elsif a_l1 < 9.504178996162932e-001 + u, v = self.pade_seven + elsif a_l1 < 2.097847961257068e+000 + u, v = self.pade_nine + else + maxnorm = 5.371920351148152 + n_squarings = {0, Math.log2(a_l1 / maxnorm).ceil.to_i}.max + a = self / 2**n_squarings + u, v = a.pade_thirteen + end + num = u + v + den = u.map(v) do |i, j| + -i + j + end + + r = den.solve(num) + + n_squarings.times do + r = r.matmul(r) + end + r + {% elsif T == Float32 %} + if a_l1 < 4.258730016922831e-001 + u, v = self.pade_three + elsif a_l1 < 1.880152677804762e+000 + u, v = self.pade_five + else + maxnorm = 3.925724783138660 + n_squarings = {0, Math.log2(a_l1 / maxnorm).ceil.to_i}.max + a = self / 2**n_squarings + u, v = a.pade_thirteen + end + num = u + v + den = u.map(v) do |i, j| + -i + j + end + + r = den.solve(num) + + n_squarings.times do + r = r.matmul(r) + end + r + {% else %} + {% raise Num::Internal::ShapeError.new("Invalid type #{T} for expm") %} + {% end %} + end + + # :nodoc: + protected def pade_three + b = [120, 60, 12, 1] + i, j = @shape + ident = Tensor(T).eye(i, j) + a2 = self.matmul(self) + + inter = a2.map(ident) do |i, j| + b[3] * i + b[1] * j + end + + u = self.matmul(inter) + v = a2.map(ident) do |i, j| + b[2] * i + b[0] * j + end + {u, v} + end + + # :nodoc: + protected def pade_five + b = [30240, 15120, 3360, 420, 30, 1] + i, j = @shape + ident = Tensor(T).eye(i, j) + a2 = self.matmul(self) + a4 = a2.matmul(a2) + + inter = a4.map(a2, ident) do |i, j, k| + b[5] * i + b[3] * j + b[1] * k + end + + u = self.matmul(inter) + v = a4.map(a2, ident) do |i, j, k| + b[4] * i + b[2] * j + b[0] * k + end + {u, v} + end + + # :nodoc: + protected def pade_seven + b = [17297280, 8648640, 1995840, 277200, 25200, 1512, 56, 1] + i, j = @shape + ident = Tensor(T).eye(i, j) + a2 = self.matmul(self) + a4 = a2.matmul(a2) + a6 = a4.matmul(a2) + + u_lhs = a6.map(a4, a2) do |i, j, k| + b[7] * i + b[5] * j + b[3] * k + end + u_rhs = ident.map { |i| b[1] * i } + + u = self.matmul(u_lhs + u_rhs) + + v_lhs = a6.map(a4, a2) do |i, j, k| + b[6] * i + b[4] * j + b[2] * k + end + + v_rhs = ident.map { |i| b[0] * i } + {u, v_lhs + v_rhs} + end + + # :nodoc: + protected def pade_nine + b = [17643225600, 8821612800, 2075673600, 302702400, 30270240, + 2162160, 110880, 3960, 90, 1] + + i, j = @shape + ident = Tensor(T).eye(i, j) + a2 = self.matmul(self) + a4 = a2.matmul(a2) + a6 = a4.matmul(a2) + a8 = a6.matmul(a2) + + u_lhs = a8.map(a6, a4) do |i, j, k| + b[9] * i + b[7] * j + b[5] * k + end + + u_rhs = a2.map(ident) do |i, j| + b[3] * i + b[1] * j + end + + u = self.matmul(u_lhs + u_rhs) + + v_lhs = a8.map(a6, a4) do |i, j, k| + b[8] * i + b[6] * j + b[4] * k + end + + v_rhs = a2.map(ident) do |i, j| + b[2] * i + b[0] * j + end + + {u, v_lhs + v_rhs} + end + + # :nodoc: + protected def pade_thirteen + b = [64764752532480000, 32382376266240000, 7771770303897600, + 1187353796428800, 129060195264000, 10559470521600, 670442572800, + 33522128640, 1323241920, 40840800, 960960, 16380, 182, 1] + + i, j = @shape + ident = Tensor(T).eye(i, j) + + a2 = self.matmul(self) + a4 = a2.matmul(a2) + a6 = a4.matmul(a2) + + u_dot_first = a6.map(a4, a2) do |i, j, k| + b[13] * i + b[11] * j + b[9] * k + end + + u_lhs = a6.map(a4, a2) do |i, j, k| + b[7] * i + b[5] * j + b[3] * k + end + + u_rhs = ident.map { |i| b[1] * i } + + u = self.matmul(a6.matmul(u_dot_first) + u_lhs + u_rhs) + + v_dot_first = a6.map(a4, a2) do |i, j, k| + b[12] * i + b[10] * j + b[8] * k + end + + v_lhs = a6.map(a4, a2) do |i, j, k| + b[6] * i + b[4] * j + b[2] * k + end + + v_rhs = ident.map { |i| b[0] * i } + + v = a6.matmul(v_dot_first) + v_lhs + v_rhs + {u, v} + end + # :nodoc: def is_matrix unless self.rank == 2 diff --git a/src/tensor/operators.cr b/src/tensor/operators.cr index 8515d5e7..c4808bf4 100644 --- a/src/tensor/operators.cr +++ b/src/tensor/operators.cr @@ -47,7 +47,7 @@ module Num end # :ditto: - def {{name}}(a : Tensor | Enumerable, b : Number) + def {{name}}(a : Tensor | Enumerable, b : Number | Complex) a_t = a.to_tensor a_t.map do |i| i {{operator.id}} b @@ -70,7 +70,7 @@ module Num end # :ditto: - def {{name}}(a : Number, b : Number) + def {{name}}(a : Number, b : Number | Complex) a {{operator.id}} b end end @@ -110,7 +110,7 @@ module Num end # :ditto: - def {{fn.id}}(a : Number) + def {{fn.id}}(a : Number | Complex) Math.{{fn.id}}(a) end end @@ -234,19 +234,19 @@ class Tensor(T) end end - def {{name}}(b : Number) + def {{name}}(b : Number | Complex) self.map do |i| i {{operator.id}} b end end - def {{operator.id}}(b : Number) + def {{operator.id}}(b : Number | Complex) self.map do |i| i {{operator.id}} b end end - def {{name}}!(b : Number) + def {{name}}!(b : Number | Complex) self.map! do |i| i {{operator.id}} b end @@ -357,4 +357,8 @@ class Tensor(T) stdlibwrap ldexp stdlibwrap max stdlibwrap min + + def - + map { |i| -i } + end end diff --git a/src/tensor/random.cr b/src/tensor/random.cr new file mode 100644 index 00000000..b5aa83fe --- /dev/null +++ b/src/tensor/random.cr @@ -0,0 +1,363 @@ +# Copyright (c) 2020 Crystal Data Contributors +# +# MIT License +# +# Permission is hereby granted, free of charge, to any person obtaining +# a copy of this software and associated documentation files (the +# "Software"), to deal in the Software without restriction, including +# without limitation the rights to use, copy, modify, merge, publish, +# distribute, sublicense, and/or sell copies of the Software, and to +# permit persons to whom the Software is furnished to do so, subject to +# the following conditions: +# +# The above copyright notice and this permission notice shall be +# included in all copies or substantial portions of the Software. +# +# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +# EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF +# MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +# NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE +# LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION +# OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION +# WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. + +class Tensor(T) + # Creates a `Tensor` sampled from a provided range, with a given + # shape. + # + # The generic types of the `Tensor` are inferred from the endpoints + # of the range + # + # Arguments + # --------- + # *r* : Range(U, U) + # Range of values to sample between + # *shape* : Array(Int) + # Shape of returned `Tensor` + # + # Examples + # -------- + # ``` + # Num::Rand.set_seed(0) + # t = Tensor.random(0...10, [2, 2]) + # t + # + # # [[8, 4], + # # [7, 4]] + # ``` + def self.random(r : Range(U, U), shape : Array(Int)) : Tensor(U) forall U + self.new(shape) do + Num::Rand.stdlib_generator.rand(r) + end + end + + # Generate random floating point values between 0 and 1 + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of Array to generate + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.rand(shape : Array(Int)) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.float + {% elsif T == Float32 %} + Num::Rand.generator.float32 + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generates a Tensor containing a beta-distribution collection + # of values + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # a : Float + # Shape parameter of distribution + # b : Float + # Shape parameter of distribution + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.beta(shape : Array(Int), a : Float, b : Float) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.beta(a: a, b: b) + {% elsif T == Float32 %} + Num::Rand.generator.beta32(a: a, b: b) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generates a Tensor containing chi-square distributed values + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # df : Float + # Degrees of freedom + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.chisq(shape : Array(Int), df : Float) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.chisq(df) + {% elsif T == Float32 %} + Num::Rand.generator.chisq32(df) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generates a Tensor containing expontentially distributed + # values + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # scale : Float = 1.0 + # Scale parameter of the distribution + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.exp(shape : Array(Int), scale : Float = 1.0) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.exp(scale) + {% elsif T == Float32 %} + Num::Rand.generator.exp32(scale) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generates a Tensor containing f-snedecor distributed + # values + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # df1 : Float + # degrees of freedom of the underlying chi-square distribution, + # numerator side; usually mentioned as m. + # df2 : Float + # degrees of freedom of the underlying chi-square distribution, + # denominator side; usually mentioned as n. + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.fsned(shape : Array(Int), df1 : Float, df2 : Float) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.f(df1, df2) + {% elsif T == Float32 %} + Num::Rand.generator.f32(df1, df2) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generate a gamma-distributed, pseudo-random Tensor + # + # Arguments + # --------- + # t_shape : Array(Int) + # Shape of output Tensor + # shape : Float + # shape parameter of the distribution; usually mentioned as k + # scale : Float = 1.0 + # scale parameter of the distribution; usually mentioned as θ + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.gamma(t_shape : Array(Int), shape : Float, scale : Float = 1.0) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.gamma(shape, scale) + {% elsif T == Float32 %} + Num::Rand.generator.gamma32(shape, scale) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generate a laplace-distributed, pseudo-random Tensor + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # loc : Float = 0.0 + # centrality parameter, or mean of the distribution; usually + # mentioned as μ + # scale : Float = 1.0 + # scale parameter of the distribution; usually mentioned as b + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.laplace(shape : Array(Int), loc : Float = 0.0, scale : Float = 1.0) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.laplace(loc, scale) + {% elsif T == Float32 %} + Num::Rand.generator.laplace32(loc, scale) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generate a log-normal-distributed, pseudo-random Tensor + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # loc : Float = 0.0 + # centrality parameter, or mean of the underlying normal distribution; + # usually mentioned as μ + # sigma : Float = 1.0 + # scale parameter, or standard deviation of the underlying normal + # distribution; usually mentioned as σ + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.lognormal(shape : Array(Int), loc : Float = 0.0, sigma : Float = 1.0) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.lognormal(loc, sigma) + {% elsif T == Float32 %} + Num::Rand.generator.lognormal32(loc, sigma) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generate a Tensor containing a normal-distribution collection + # of values + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of Tensor to create + # loc = 0.0 + # Centrality parameter + # sigma = 1.0 + # Standard deviation + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.normal(shape : Array(Int), loc = 0.0, sigma = 1.0) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.normal(loc, sigma) + {% elsif T == Float32 %} + Num::Rand.generator.normal32(loc, sigma) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end + + # Generate a poisson-distributed, pseudo-random Tensor(Int64) + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # lam : Float = 1.0 + # separation parameter of the distribution; usually mentioned as λ + # + # Returns + # ------- + # Tensor(Int64) + # + # Examples + # -------- + def self.poisson(shape : Array(Int), lam : Float = 1.0) : Tensor(Int64) + self.new(shape) do + Num::Rand.generator.poisson(lam) + end + end + + # Generate a t-student-distributed, pseudo-random Tensor + # + # Arguments + # --------- + # shape : Array(Int) + # Shape of output Tensor + # df : Float + # degrees of freedom of the distribution; usually mentioned as n + # + # Returns + # ------- + # Tensor(T) + # + # Examples + # -------- + def self.t_student(shape : Array(Int), df : Float) : Tensor(T) + self.new(shape) do + {% if T == Float64 %} + Num::Rand.generator.t(df) + {% elsif T == Float32 %} + Num::Rand.generator.t32(df) + {% else %} + {% raise "Invalid dtype #{T} for random methods" %} + {% end %} + end + end +end diff --git a/src/tensor/reductions.cr b/src/tensor/reductions.cr index c6618911..d8a8bf91 100644 --- a/src/tensor/reductions.cr +++ b/src/tensor/reductions.cr @@ -488,6 +488,64 @@ module Num end end + # Sorts a `Tensor`, treating it's elements like the `Tensor` + # is flat. + # + # Arguments + # --------- + # *a* : Tensor | Enumerable + # Argument to sort + # *axis* : Int + # + # Examples + # -------- + # ``` + # a = [3, 2, 1].to_tensor + # Num.sort(a) # => [1, 2, 3] + # ``` + def sort(a : Tensor | Enumerable) + a_t = a.to_tensor + ret = a_t.dup(Num::RowMajor) + Slice.new(ret.to_unsafe, ret.size).sort! + ret + end + + # Sorts a `Tensor` along an axis. + # + # Arguments + # --------- + # *a* : Tensor | Enumerable + # Argument to sort + # *axis* : Int + # Axis to sort along + # + # Examples + # -------- + # ``` + # t = Tensor.random(0...10, [3, 3, 2]) + # puts Num.sort(t, axis: 1) + # + # # [[[1, 1], + # # [4, 5], + # # [5, 7]], + # # + # # [[0, 0], + # # [2, 3], + # # [8, 4]], + # # + # # [[2, 5], + # # [5, 7], + # # [5, 7]]] + # ``` + def sort(a : Tensor | Enumerable, axis : Int) + a_t = a.to_tensor + ret = a_t.dup(Num::RowMajor) + ret.yield_along_axis(axis) do |view| + view[...] = Num.sort(view) + end + ret + end + # Asserts that two `Tensor`s are equal, allowing for small # margins of errors with floating point values using # an EPSILON value. @@ -516,7 +574,9 @@ module Num unless a.size == b.size return false end - a.map(b) do |i, j| + a_t = a.to_tensor + b_t = b.to_tensor + a_t.map(b_t) do |i, j| m = (i - j).abs < epsilon unless m return false @@ -580,6 +640,22 @@ module Num ) end + # Returns a hash containing the count of each + # unique element of a `Tensor` + # + # Arguments + # --------- + # *a* : Tensor + # Tensor to count + # *axis* : Tens + # + # + # Examples + # -------- + # ``` + # a = [[3, 4], [2, 2]] + # Num.value_counts(a) # => {3 => 1, 4 => 1, 2 => 2} + # ``` def value_counts(a : Tensor(U)) forall U counts = Hash(U, Int32).new a.each do |e| @@ -1024,6 +1100,72 @@ class Tensor(T) Num.ptp(self, axis, dims) end + # Sorts a `Tensor`, treating it's elements like the `Tensor` + # is flat. + # + # Arguments + # --------- + # *a* : Tensor | Enumerable + # Argument to sort + # *axis* : Int + # + # Examples + # -------- + # ``` + # a = [3, 2, 1].to_tensor + # Num.sort(a) # => [1, 2, 3] + # ``` + def sort + Num.sort(self) + end + + # Sorts a `Tensor` along an axis. + # + # Arguments + # --------- + # *a* : Tensor | Enumerable + # Argument to sort + # *axis* : Int + # Axis to sort along + # + # Examples + # -------- + # ``` + # t = Tensor.random(0...10, [3, 3, 2]) + # puts Num.sort(t, axis: 1) + # + # # [[[1, 1], + # # [4, 5], + # # [5, 7]], + # # + # # [[0, 0], + # # [2, 3], + # # [8, 4]], + # # + # # [[2, 5], + # # [5, 7], + # # [5, 7]]] + # ``` + def sort(axis : Int) + Num.sort(self, axis) + end + + # Returns a hash containing the count of each + # unique element of a `Tensor` + # + # Arguments + # --------- + # *a* : Tensor + # Tensor to count + # *axis* : Tens + # + # + # Examples + # -------- + # ``` + # a = [[3, 4], [2, 2]] + # Num.value_counts(a) # => {3 => 1, 4 => 1, 2 => 2} + # ``` def value_counts Num.value_counts(self) end diff --git a/src/tensor/tensor.cr b/src/tensor/tensor.cr index 71d91b13..ffc61517 100644 --- a/src/tensor/tensor.cr +++ b/src/tensor/tensor.cr @@ -26,6 +26,7 @@ require "./internal/constants" require "./internal/utils" require "./internal/print" require "./internal/iteration" +require "./internal/yield_iterators" require "./internal/ndindex" require "../libs/cblas" @@ -467,6 +468,42 @@ class Tensor(T) slice(args) end + # Advanced indexing operation for a `Tensor`, allows Tensor + # to be sliced in a non-viewable manner, creating a copy + # but allowing greater flexibility. + # + # For now, only arrays of integers are supported. + # + # Arguments + # --------- + # args : Array(Int) + # Array containing the indexers in each dimension + # + # Returns + # ------- + # Tensor(T) + # A copy of the sliced data + # + # Examples + # -------- + # ``` + # t = Tensor.new([3, 3, 2]) { |i| i } + # t[[0, 2, 1], [1, 0, 1]]? + # + # # [[0, 1], + # # [2, 3], + # # [6, 7]] + # ``` + # + def []?(*args : Array(Int)) + copy_slice(args.to_a) + end + + # :ditto: + def []?(args : Array(Array(Int))) + copy_slice(args) + end + # The primary method of setting Tensor values. The slicing behavior # for this method is identical to the `[]` method. # @@ -615,7 +652,7 @@ class Tensor(T) # # 3 # ``` def each - iter.each do |el| + strided_iteration(self) do |_, el| yield el.value end end @@ -644,7 +681,7 @@ class Tensor(T) # # 3 # ``` def each_pointer - iter.each do |el| + strided_iteration(self) do |_, el| yield el end end @@ -809,8 +846,9 @@ class Tensor(T) # ``` def map(&block : T -> U) : Tensor(U) forall U t = Tensor(U).new(@shape) - t.iter(self).each do |i, j| - i.value = yield j.value + data = t.to_unsafe + strided_iteration(self) do |index, el| + data[index] = yield el.value end t end @@ -862,8 +900,9 @@ class Tensor(T) def map(t : Tensor(U), &block : T, U -> V) : Tensor(V) forall U, V a, b = Num::Internal.broadcast(self, t) r = Tensor(V).new(a.shape) - r.iter(a, b).each do |i, j, k| - i.value = yield(j.value, k.value) + data = r.to_unsafe + dual_strided_iteration(a, b) do |index, a, b| + data[index] = yield a.value, b.value end r end @@ -895,7 +934,7 @@ class Tensor(T) # ``` def map!(t : Tensor, &block) t = t.as_shape(@shape) - iter(t).each do |i, j| + dual_strided_iteration(self, t) do |_, i, j| type_inference i, i, j end end @@ -933,8 +972,9 @@ class Tensor(T) ) : Tensor(W) forall U, V, W a, b, c = Num::Internal.broadcast(self, t, v) r = Tensor(W).new(a.shape) - r.iter(a, b, c).each do |i, j, k, l| - i.value = yield(j.value, k.value, l.value) + data = r.to_unsafe + tri_strided_iteration(a, b, c) do |index, i, j, k| + data[index] = yield i.value, j.value, k.value end r end @@ -968,10 +1008,10 @@ class Tensor(T) # a.map!(b, c) { |i, j, k| i + j + k } # a # => [0, 3, 6] # ``` - def map!(t : Tensor, v : Tensor) + def map!(t : Tensor, v : Tensor, &block) t = t.as_shape(@shape) v = v.as_shape(@shape) - iter(t, v).each do |i, j, k| + tri_strided_iteration(self, t, v) do |index, i, j, k| type_inference i, i, j, k end end @@ -1118,7 +1158,11 @@ class Tensor(T) def as_type(u : U.class) : Tensor(U) forall U r = Tensor(U).new(@shape) r.map!(self) do |_, j| - j + {% if T == Bool %} + j ? 1 : 0 + {% else %} + j + {% end %} end r end @@ -1429,6 +1473,10 @@ class Tensor(T) end end + def iter_attrs + {@shape.to_unsafe, @strides.to_unsafe, self.rank} + end + # :nodoc: def map_along_axis(axis : Int) if axis < 0 @@ -1618,6 +1666,28 @@ class Tensor(T) true end + private def copy_slice(args : Array(Array(Int))) + s0 = args[0].size + args.each do |arg| + unless s0 == arg.size + raise Num::Internal::ShapeError.new("All indexers must be the same size") + end + end + + unless s0 <= self.rank + raise Num::Internal::ShapeError.new("Too many dimensions for Tensor") + end + + out_shape = [s0] + @shape[args.size...] + + ret = Tensor(T).new(out_shape) + s0.times do |i| + idx = args.map { |arg| arg[i] } + ret[i] = self[idx] + end + ret + end + private def slice(args : Array) new_shape = @shape.dup new_strides = @strides.dup @@ -1670,6 +1740,12 @@ class Tensor(T) end private def normalize(arg : Range, i : Int32) + a_end = arg.end + if a_end.is_a?(Int32) + if a_end > @shape[i] + arg = arg.begin...@shape[i] + end + end s, o = Indexable.range_to_index_and_count(arg, @shape[i]) if s >= @shape[i] raise Num::Internal::IndexError.new(