Unverified 提交 f772066a authored 作者: Jesse Grabowski's avatar Jesse Grabowski 提交者: GitHub

Remove `sparse.sandbox` (#1664)

* Delete `sparse.sandbox` * delete test_sp2.py
上级 95acdb3d
"""
Convolution-like operations with sparse matrix multiplication.
To read about different sparse formats, see
U{http://www-users.cs.umn.edu/~saad/software/SPARSKIT/paper.ps}.
@todo: Automatic methods for determining best sparse format?
"""
# COPIED FROM hpu/icml09/sp.py
import numpy as np
from scipy import sparse as scipy_sparse
import pytensor
import pytensor.sparse
from pytensor import sparse
from pytensor import tensor as pt
from pytensor.graph.op import Op
from pytensor.tensor.math import dot
from pytensor.tensor.math import max as pt_max
from pytensor.tensor.shape import reshape
def register_specialize(lopt, *tags, **kwargs):
pytensor.compile.optdb["specialize"].register(
(kwargs and kwargs.pop("name")) or lopt.__name__, lopt, "fast_run", *tags
)
class ConvolutionIndices(Op):
"""Build indices for a sparse CSC matrix that could implement A
(convolve) B.
This generates a sparse matrix M, which generates a stack of
image patches when computing the dot product of M with image
patch. Convolution is then simply the dot product of (img x M)
and the kernels.
"""
__props__ = ()
@staticmethod
def conv_eval(inshp, kshp, strides=(1, 1), mode="valid"):
(dx, dy) = strides
return convolution_indices.evaluate(inshp, kshp, (dx, dy), mode=mode, ws=True)
# img_shape and ker_shape are (height,width)
@staticmethod
def evaluate(inshp, kshp, strides=(1, 1), nkern=1, mode="valid", ws=True):
"""Build a sparse matrix which can be used for performing...
* convolution: in this case, the dot product of this matrix
with the input images will generate a stack of images
patches. Convolution is then a tensordot operation of the
filters and the patch stack.
* sparse local connections: in this case, the sparse matrix
allows us to operate the weight matrix as if it were
fully-connected. The structured-dot with the input image gives
the output for the following layer.
:param ker_shape: shape of kernel to apply (smaller than image)
:param img_shape: shape of input images
:param mode: 'valid' generates output only when kernel and
image overlap overlap fully. Convolution obtained
by zero-padding the input
:param ws: must be always True
:param (dx,dy): offset parameter. In the case of no weight sharing,
gives the pixel offset between two receptive fields.
With weight sharing gives the offset between the
top-left pixels of the generated patches
:rtype: tuple(indices, indptr, logical_shape, sp_type, out_img_shp)
:returns: the structure of a sparse matrix, and the logical dimensions
of the image which will be the result of filtering.
"""
if not ws:
raise Exception("ws is obsolete and it must be always True")
(dx, dy) = strides
# inshp contains either 2 entries (height,width) or 3 (nfeatures,h,w)
# in the first case, default nfeatures to 1
if np.size(inshp) == 2:
inshp = (1, *inshp)
inshp = np.array(inshp)
kshp = np.array(kshp)
ksize = np.prod(kshp)
# kern = ksize - 1 - np.arange(ksize)
# size of output image if doing proper convolution
# (mode='full',dx=dy=0) outshp is the actual output shape
# given the parameters
fulloutshp = inshp[1:] + kshp - 1
if mode == "valid":
s = -1
else:
s = 1
outshp = np.int64(
np.ceil((inshp[1:] + s * kshp - s * 1) / np.array([dy, dx], dtype="float"))
)
if any(outshp <= 0):
err = (
"Invalid kernel",
kshp,
"and/or step size",
(dx, dy),
"for given input shape",
inshp,
)
raise ValueError(err)
outsize = np.prod(outshp)
insize = np.prod(inshp)
# range of output units over which to iterate
if mode == "valid":
lbound = np.array([kshp[0] - 1, kshp[1] - 1])
ubound = lbound + (inshp[1:] - kshp + 1)
else:
lbound = np.zeros(2)
ubound = fulloutshp
# coordinates of image in "fulloutshp" coordinates
topleft = np.array([kshp[0] - 1, kshp[1] - 1])
# bound when counting the receptive field
botright = topleft + inshp[1:]
# sparse matrix specifics...
if ws:
spmatshp = (outsize * np.prod(kshp) * inshp[0], insize)
else:
spmatshp = (nkern * outsize, insize)
spmat = scipy_sparse.lil_matrix(spmatshp)
# loop over output image pixels
# z, zz = 0, 0
# incremented every time we write something to the sparse
# matrix this is used to track the ordering of filter tap
# coefficient in sparse column ordering
tapi, ntaps = 0, 0
# Note: looping over the number of kernels could've been done
# more efficiently as the last step (when writing to
# spmat). However, this messes up the ordering of the column
# values (order in which you write the values determines how
# the vectorized data will get used later one)
for fmapi in range(inshp[0]): # loop over input features
# loop over number of kernels (nkern=1 for weight sharing)
for n in range(nkern):
# FOR EACH OUTPUT PIXEL...
# loop over output image height
for oy in np.arange(lbound[0], ubound[0], dy, dtype=int):
# loop over output image width
for ox in np.arange(lbound[1], ubound[1], dx, dtype=int):
# kern[l] is filter value to apply at (oj,oi)
# for (iy,ix)
l = 0
# ... ITERATE OVER INPUT UNITS IN RECEPTIVE FIELD
for ky in oy + np.arange(kshp[0], dtype=int):
for kx in ox + np.arange(kshp[1], dtype=int):
# verify if we are still within image
# boundaries. Equivalent to
# zero-padding of the input image
if all((ky, kx) >= topleft) and all(
(ky, kx) < botright
):
# convert to "valid" input space
# coords used to determine column
# index to write to in sparse mat
iy, ix = np.array((ky, kx), dtype=int) - topleft
# determine raster-index of input pixel...
# taking into account multiple
# input features
col = int(
iy * inshp[2]
+ ix
+ fmapi * np.prod(inshp[1:], dtype=int)
)
# convert oy,ox values to output
# space coordinates
if mode == "full":
(y, x) = (oy, ox)
else:
(y, x) = (oy, ox) - topleft
# taking into account step size
(y, x) = np.array([y, x], dtype=int) / (dy, dx)
# convert to row index of sparse matrix
if ws:
row = int(
(y * outshp[1] + x) * inshp[0] * ksize
+ l
+ fmapi * ksize
)
else:
row = int(y * outshp[1] + x)
# Store something at that location
# in sparse matrix. The written
# value is only useful for the
# sparse case. It will determine
# the way kernel taps are mapped
# onto the sparse columns (idea of
# kernel map)
# n*... only for sparse
spmat[int(row + n * outsize), int(col)] = tapi + 1
# total number of active taps
# (used for kmap)
ntaps += 1
# absolute tap index (total number of taps)
tapi += 1
# move on to next filter tap l=(l+1)%ksize
l += 1
if spmat.format != "csc":
spmat = spmat.tocsc().sorted_indices()
else:
# BUG ALERT: scipy0.6 has bug where data and indices are written in
# reverse column ordering.
# Explicit call to sorted_indices removes this problem.
spmat = spmat.sorted_indices()
if ws:
kmap = None
else:
kmap = np.zeros(ntaps, dtype="int")
k = 0
# print 'TEMPORARY BUGFIX: REMOVE !!!'
for j in range(spmat.shape[1]):
for i_idx in range(spmat.indptr[j], spmat.indptr[j + 1]):
if spmat.data[i_idx] != 0:
# this is == spmat[i,j] - 1
kmap[k] = spmat.data[i_idx] - 1
k += 1
# when in valid mode, it is more efficient to store in sparse row
# TODO: need to implement structured dot for csr matrix
assert spmat.format == "csc"
sptype = "csc"
# sptype = 'csr' if mode=='valid' else 'csc'
if 0 and mode == "valid":
spmat = spmat.tocsr()
rval = (spmat.indices[: spmat.size], spmat.indptr, spmatshp, sptype, outshp)
if kmap is not None:
rval += (kmap,)
return rval
def perform(self, node, inputs, outputs):
(inshp, kshp) = inputs
(out_indices, out_indptr, spmat_shape) = outputs
indices, indptr, spmatshp, _outshp = self.evaluate(inshp, kshp)
out_indices[0] = indices
out_indptr[0] = indptr
spmat_shape[0] = np.asarray(spmatshp)
convolution_indices = ConvolutionIndices()
def convolve(
kerns,
kshp,
nkern,
images,
imgshp,
step=(1, 1),
bias=None,
mode="valid",
flatten=True,
):
"""Convolution implementation by sparse matrix multiplication.
:note: For best speed, put the matrix which you expect to be
smaller as the 'kernel' argument
"images" is assumed to be a matrix of shape batch_size x img_size,
where the second dimension represents each image in raster order
If flatten is "False", the output feature map will have shape:
.. code-block:: python
batch_size x number of kernels x output_size
If flatten is "True", the output feature map will have shape:
.. code-block:: python
batch_size x number of kernels * output_size
.. note::
IMPORTANT: note that this means that each feature map (image
generate by each kernel) is contiguous in memory. The memory
layout will therefore be: [ <feature_map_0> <feature_map_1>
... <feature_map_n>], where <feature_map> represents a
"feature map" in raster order
kerns is a 2D tensor of shape nkern x N.prod(kshp)
:param kerns: 2D tensor containing kernels which are applied at every pixel
:param kshp: tuple containing actual dimensions of kernel (not symbolic)
:param nkern: number of kernels/filters to apply.
nkern=1 will apply one common filter to all input pixels
:param images: tensor containing images on which to apply convolution
:param imgshp: tuple containing image dimensions
:param step: determines number of pixels between adjacent receptive fields
(tuple containing dx,dy values)
:param mode: 'full', 'valid' see CSM.evaluate function for details
:param sumdims: dimensions over which to sum for the tensordot operation.
By default ((2,),(1,)) assumes kerns is a nkern x kernsize
matrix and images is a batchsize x imgsize matrix
containing flattened images in raster order
:param flatten: flatten the last 2 dimensions of the output. By default,
instead of generating a batchsize x outsize x nkern tensor,
will flatten to batchsize x outsize*nkern
:return: out1, symbolic result
:return: out2, logical shape of the output img (nkern,height,width)
:TODO: test for 1D and think of how to do n-d convolutions
"""
# start by computing output dimensions, size, etc
kern_size = np.int64(np.prod(kshp))
# inshp contains either 2 entries (height,width) or 3 (nfeatures,h,w)
# in the first case, default nfeatures to 1
if np.size(imgshp) == 2:
imgshp = (1, *imgshp)
# construct indices and index pointers for sparse matrix, which,
# when multiplied with input images will generate a stack of image
# patches
indices, indptr, spmat_shape, sptype, outshp = convolution_indices.conv_eval(
imgshp, kshp, step, mode
)
# build sparse matrix, then generate stack of image patches
csc = pytensor.sparse.CSM(sptype)(
np.ones(indices.size), indices, indptr, spmat_shape
)
patches = (sparse.structured_dot(csc, images.T)).T
# compute output of linear classifier
pshape = pt.stack(
[
images.shape[0] * pt.as_tensor(np.prod(outshp)),
pt.as_tensor(imgshp[0] * kern_size),
]
)
patch_stack = reshape(patches, pshape, ndim=2)
# kern is of shape: nkern x ksize*number_of_input_features
# output is thus of shape: bsize*outshp x nkern
output = dot(patch_stack, kerns.T)
# add bias across each feature map (more efficient to do it now)
if bias is not None:
output += bias
# now to have feature maps in raster order ...
# go from bsize*outshp x nkern to bsize x nkern*outshp
newshp = pt.stack(
[images.shape[0], pt.as_tensor(np.prod(outshp)), pt.as_tensor(nkern)]
)
tensout = reshape(output, newshp, ndim=3)
output = tensout.transpose(0, 2, 1)
if flatten:
output = pt.flatten(output, 2)
return output, np.hstack((nkern, outshp))
def max_pool(images, imgshp, maxpoolshp):
"""Implements a max pooling layer
Takes as input a 2D tensor of shape batch_size x img_size and
performs max pooling. Max pooling downsamples by taking the max
value in a given area, here defined by maxpoolshp. Outputs a 2D
tensor of shape batch_size x output_size.
:param images: 2D tensor containing images on which to apply convolution.
Assumed to be of shape batch_size x img_size
:param imgshp: tuple containing image dimensions
:param maxpoolshp: tuple containing shape of area to max pool over
:return: out1, symbolic result (2D tensor)
:return: out2, logical shape of the output
"""
poolsize = np.int64(np.prod(maxpoolshp))
# imgshp contains either 2 entries (height,width) or 3 (nfeatures,h,w)
# in the first case, default nfeatures to 1
if np.size(imgshp) == 2:
imgshp = (1, *imgshp)
# construct indices and index pointers for sparse matrix, which,
# when multiplied with input images will generate a stack of image
# patches
indices, indptr, spmat_shape, sptype, outshp = convolution_indices.conv_eval(
imgshp, maxpoolshp, maxpoolshp, mode="valid"
)
# print 'XXXXXXXXXXXXXXXX MAX POOLING LAYER XXXXXXXXXXXXXXXXXXXX'
# print 'imgshp = ', imgshp
# print 'maxpoolshp = ', maxpoolshp
# print 'outshp = ', outshp
# build sparse matrix, then generate stack of image patches
csc = pytensor.sparse.CSM(sptype)(
np.ones(indices.size), indices, indptr, spmat_shape
)
patches = sparse.structured_dot(csc, images.T).T
pshape = pt.stack(
[
images.shape[0] * pt.as_tensor(np.prod(outshp)),
pt.as_tensor(imgshp[0]),
pt.as_tensor(poolsize),
]
)
patch_stack = reshape(patches, pshape, ndim=3)
out1 = pt_max(patch_stack, axis=2)
pshape = pt.stack(
[
images.shape[0],
pt.as_tensor(np.prod(outshp)),
pt.as_tensor(imgshp[0]),
]
)
out2 = reshape(out1, pshape, ndim=3)
out3 = out2.transpose(0, 2, 1)
return pt.flatten(out3, 2), outshp
import numpy as np
import scipy.sparse
import pytensor
from pytensor import tensor as pt
from pytensor.graph.basic import Apply
from pytensor.graph.op import Op
from pytensor.sparse.basic import (
Remove0,
SparseTensorType,
_is_sparse,
as_sparse_variable,
remove0,
)
from pytensor.tensor.type import discrete_dtypes, float_dtypes
# Probability Ops are currently back in sandbox, because they do not respect
# PyTensor's Op contract, as their behaviour is not reproducible: calling
# the perform() method twice with the same argument will yield different
# results.
# from pytensor.sparse.basic import (
# Multinomial, multinomial, Poisson, poisson,
# Binomial, csr_fbinomial, csc_fbinomial, csr_dbinomial, csc_dbinomial)
# Alias to maintain compatibility
EliminateZeros = Remove0
eliminate_zeros = remove0
# Probability
class Poisson(Op):
"""Return a sparse having random values from a Poisson density
with mean from the input.
WARNING: This Op is NOT deterministic, as calling it twice with the
same inputs will NOT give the same result. This is a violation of
PyTensor's contract for Ops
:param x: Sparse matrix.
:return: A sparse matrix of random integers of a Poisson density
with mean of `x` element wise.
"""
__props__ = ()
def make_node(self, x):
x = as_sparse_variable(x)
return Apply(self, [x], [x.type()])
def perform(self, node, inputs, outputs):
(x,) = inputs
(out,) = outputs
assert _is_sparse(x)
assert x.format in ("csr", "csc")
out[0] = x.copy()
out[0].data = np.asarray(np.random.poisson(out[0].data), dtype=x.dtype)
out[0].eliminate_zeros()
def grad(self, inputs, outputs_gradients):
comment = "No gradient exists for class Poisson in\
pytensor/sparse/sandbox/sp2.py"
return [
pytensor.gradient.grad_undefined(
op=self, x_pos=0, x=inputs[0], comment=comment
)
]
def infer_shape(self, fgraph, node, ins_shapes):
return ins_shapes
poisson = Poisson()
class Binomial(Op):
"""Return a sparse matrix having random values from a binomial
density having number of experiment `n` and probability of success
`p`.
WARNING: This Op is NOT deterministic, as calling it twice with the
same inputs will NOT give the same result. This is a violation of
PyTensor's contract for Ops
:param n: Tensor scalar representing the number of experiment.
:param p: Tensor scalar representing the probability of success.
:param shape: Tensor vector for the output shape.
:return: A sparse matrix of integers representing the number
of success.
"""
__props__ = ("format", "dtype")
def __init__(self, format, dtype):
self.format = format
self.dtype = np.dtype(dtype).name
def make_node(self, n, p, shape):
n = pt.as_tensor_variable(n)
p = pt.as_tensor_variable(p)
shape = pt.as_tensor_variable(shape)
assert n.dtype in discrete_dtypes
assert p.dtype in float_dtypes
assert shape.dtype in discrete_dtypes
return Apply(
self,
[n, p, shape],
[SparseTensorType(dtype=self.dtype, format=self.format)()],
)
def perform(self, node, inputs, outputs):
(n, p, shape) = inputs
(out,) = outputs
binomial = np.random.binomial(n, p, size=shape)
csx_matrix = getattr(scipy.sparse, self.format + "_matrix")
out[0] = csx_matrix(binomial, dtype=self.dtype)
def connection_pattern(self, node):
return [[True], [True], [False]]
def grad(self, inputs, gout):
(n, p, _shape) = inputs
(_gz,) = gout
comment_n = "No gradient exists for the number of samples in class\
Binomial of pytensor/sparse/sandbox/sp2.py"
comment_p = "No gradient exists for the prob of success in class\
Binomial of pytensor/sparse/sandbox/sp2.py"
return [
pytensor.gradient.grad_undefined(op=self, x_pos=0, x=n, comment=comment_n),
pytensor.gradient.grad_undefined(op=self, x_pos=1, x=p, comment=comment_p),
pytensor.gradient.disconnected_type(),
]
def infer_shape(self, fgraph, node, ins_shapes):
return [(node.inputs[2][0], node.inputs[2][1])]
csr_fbinomial = Binomial("csr", "float32")
csc_fbinomial = Binomial("csc", "float32")
csr_dbinomial = Binomial("csr", "float64")
csc_dbinomial = Binomial("csc", "float64")
class Multinomial(Op):
"""Return a sparse matrix having random values from a multinomial
density having number of experiment `n` and probability of success
`p`.
WARNING: This Op is NOT deterministic, as calling it twice with the
same inputs will NOT give the same result. This is a violation of
PyTensor's contract for Ops
:param n: Tensor type vector or scalar representing the number of
experiment for each row. If `n` is a scalar, it will be
used for each row.
:param p: Sparse matrix of probability where each row is a probability
vector representing the probability of success. N.B. Each row
must sum to one.
:return: A sparse matrix of random integers from a multinomial density
for each row.
:note: It will works only if `p` have csr format.
"""
__props__ = ()
def make_node(self, n, p):
n = pt.as_tensor_variable(n)
p = as_sparse_variable(p)
assert p.format in ("csr", "csc")
return Apply(self, [n, p], [p.type()])
def perform(self, node, inputs, outputs):
(n, p) = inputs
(out,) = outputs
assert _is_sparse(p)
if p.format != "csr":
raise NotImplementedError
out[0] = p.copy()
if n.ndim == 0:
for i in range(p.shape[0]):
k, l = p.indptr[i], p.indptr[i + 1]
out[0].data[k:l] = np.random.multinomial(n, p.data[k:l])
elif n.ndim == 1:
if n.shape[0] != p.shape[0]:
raise ValueError(
"The number of element of n must be "
"the same as the number of row of p."
)
for i in range(p.shape[0]):
k, l = p.indptr[i], p.indptr[i + 1]
out[0].data[k:l] = np.random.multinomial(n[i], p.data[k:l])
def grad(self, inputs, outputs_gradients):
comment_n = "No gradient exists for the number of samples in class\
Multinomial of pytensor/sparse/sandbox/sp2.py"
comment_p = "No gradient exists for the prob of success in class\
Multinomial of pytensor/sparse/sandbox/sp2.py"
return [
pytensor.gradient.grad_undefined(
op=self, x_pos=0, x=inputs[0], comment=comment_n
),
pytensor.gradient.grad_undefined(
op=self, x_pos=1, x=inputs[1], comment=comment_p
),
]
def infer_shape(self, fgraph, node, ins_shapes):
return [ins_shapes[1]]
multinomial = Multinomial()
import time
import pytest
pytest.importorskip("scipy", minversion="0.7.0")
import numpy as np
from scipy.signal import convolve2d
from pytensor import function
from pytensor.sparse.sandbox import sp
from pytensor.tensor.type import dmatrix, dvector
from tests import unittest_tools as utt
class TestSP:
@pytest.mark.slow
def test_convolution(self):
# print '\n\n*************************************************'
# print ' TEST CONVOLUTION'
# print '*************************************************'
# fixed parameters
bsize = 10 # batch size
imshp = (28, 28)
kshp = (5, 5)
nkern = 5
ssizes = ((1, 1), (2, 2), (3, 3), (4, 4))
convmodes = ("full", "valid")
# symbolic stuff
bias = dvector()
kerns = dmatrix()
input = dmatrix()
rng = np.random.default_rng(3423489)
filters = rng.standard_normal((nkern, np.prod(kshp)))
biasvals = rng.standard_normal(nkern)
for mode in ("FAST_COMPILE", "FAST_RUN"):
ttot, ntot = 0, 0
for conv_mode in convmodes:
for ss in ssizes:
output, _outshp = sp.convolve(
kerns, kshp, nkern, input, imshp, ss, bias=bias, mode=conv_mode
)
f = function([kerns, bias, input], output, mode=mode)
# now test with real values
img2d = np.arange(bsize * np.prod(imshp)).reshape((bsize, *imshp))
img1d = img2d.reshape(bsize, -1)
# create filters (need to be flipped to use convolve2d)
filtersflipped = np.zeros((nkern, *kshp))
for k in range(nkern):
it = reversed(filters[k, :])
for i in range(kshp[0]):
for j in range(kshp[1]):
filtersflipped[k, i, j] = next(it)
# compute output with convolve2d
if conv_mode == "valid":
fulloutshp = np.array(imshp) - np.array(kshp) + 1
else:
fulloutshp = np.array(imshp) + np.array(kshp) - 1
ntime1 = time.perf_counter()
refout = np.zeros((bsize, *fulloutshp, nkern))
for b in range(bsize):
for n in range(nkern):
refout[b, ..., n] = convolve2d(
img2d[b, :, :], filtersflipped[n, ...], conv_mode
)
ntot += time.perf_counter() - ntime1
# need to flatten images
bench1 = refout[:, 0 :: ss[0], 0 :: ss[1], :].reshape(
bsize, -1, nkern
)
bench1 += biasvals.reshape(1, 1, nkern)
# swap the last two dimensions (output needs to be nkern x outshp)
bench1 = np.swapaxes(bench1, 1, 2)
ttime1 = time.perf_counter()
out1 = f(filters, biasvals, img1d)
ttot += time.perf_counter() - ttime1
temp = bench1.flatten() - out1.flatten()
assert (temp < 1e-5).all()
# test downward propagation -- symbolic stuff
# vis = pytensor.gradient.grad(output, input, output)
# downprop = function([kerns,input], vis, mode=mode)
# visval = downprop(filters,img1d)
# test downward propagation -- reference implementation
# pshape = (img1d.shape[0],np.prod(outshp[1:]),np.prod(kshp))
# patchstack = np.zeros(pshape)
# for bi in np.arange(pshape[0]): # batch index
# abspos = 0
# for outy in np.arange(outshp[1]):
# for outx in np.arange(outshp[2]):
# for ni in np.arange(nkern):
# print 'filters[n,:].shape = ', filters[n,:].shape
# print 'out1[bi,abspos].shape =',out1[bi,abspos].shape
# patchstack[bi,abspos,:] = filters[n,:]*out1[bi,abspos]
# abspos+=1
# patchstack = patchstack.reshape(1,-1)
# indices, indptr, spmat_shape, sptype, outshp = \
# sp.convolution_indices.conv_eval(imshp,kshp,ss,conv_mode)
# spmat = sparse.csc_matrix((np.ones_like(indices),indices,indptr),spmat_shape)
# visref = np.dot(patchstack, spmat.todense())
# print 'visval = ', visval
# print 'visref = ', visref
# assert np.all(visref==visval)
# print '**** Convolution Profiling Results (',mode,') ****'
# print 'Numpy processing time: ', ntot
# print 'PyTensor processing time: ', ttot
# this doesn't compare the output of anything... but I manually verified that the patches
# are properly generated
def test_multilayer_conv(self):
# fixed parameters
bsize = 10 # batch size
imshp = (5, 5)
kshp = ((3, 3), (2, 2))
nkerns = (3, 6) # per output pixel
ssizes = (((1, 1), (2, 2)),)
convmodes = ("full",) # 'valid',)
# symbolic stuff
kerns = [dmatrix(), dmatrix()]
input = dmatrix()
# build actual input images
img2d = np.arange(bsize * np.prod(imshp)).reshape((bsize, *imshp))
img1d = img2d.reshape(bsize, -1)
for mode in ("FAST_COMPILE", "FAST_RUN"):
for conv_mode in convmodes:
for ss in ssizes:
l1hid, l1shp = sp.convolve(
kerns[0],
kshp[0],
nkerns[0],
input,
imshp,
ss[0],
mode=conv_mode,
)
l1propup = function([kerns[0], input], l1hid, mode=mode)
l1kernvals = np.arange(nkerns[0] * np.prod(kshp[0])).reshape(
nkerns[0], np.prod(kshp[0])
)
l1hidval = l1propup(l1kernvals, img1d)
# actual values
l2hid, _l2shp = sp.convolve(
kerns[1],
kshp[1],
nkerns[1],
l1hid,
l1shp,
ss[1],
mode=conv_mode,
)
l2propup = function([kerns[1], l1hid], l2hid, mode=mode)
l2kernvals = np.arange(
nkerns[1] * np.prod(kshp[1]) * nkerns[0]
).reshape(nkerns[1], np.prod(kshp[1]) * nkerns[0])
# for debugging, we bring things back to integers
l1hidval = np.arange(np.size(l1hidval)).reshape(l1hidval.shape)
l2propup(l2kernvals, l1hidval)
def test_maxpool(self):
# generate flatted images
maxpoolshps = ((2, 2), (3, 3), (4, 4), (5, 5), (6, 6))
rng = np.random.default_rng(2938)
imval = rng.random((4, 5, 10, 10))
images = dmatrix()
for maxpoolshp in maxpoolshps:
# symbolic stuff
output, outshp = sp.max_pool(images, imval.shape[1:], maxpoolshp)
f = function(
[
images,
],
[
output,
],
)
output_val = f(imval.reshape(imval.shape[0], -1))
# numeric verification
my_output_val = np.zeros(
(
imval.shape[0],
imval.shape[1],
imval.shape[2] // maxpoolshp[0],
imval.shape[3] // maxpoolshp[1],
)
)
assert np.prod(my_output_val.shape[1:]) == np.prod(
np.r_[imval.shape[1], outshp]
)
for n in range(imval.shape[0]):
for k in range(imval.shape[1]):
for i in range(imval.shape[2] // maxpoolshp[0]):
for j in range(imval.shape[3] // maxpoolshp[1]):
ii, jj = i * maxpoolshp[0], j * maxpoolshp[1]
patch = imval[
n, k, ii : ii + maxpoolshp[0], jj : jj + maxpoolshp[1]
]
my_output_val[n, k, i, j] = np.max(patch)
my_output_val = my_output_val.reshape(imval.shape[0], -1)
assert np.all(output_val == my_output_val)
def mp(input):
output, _outshp = sp.max_pool(input, imval.shape[1:], maxpoolshp)
return output
utt.verify_grad(mp, [imval.reshape(imval.shape[0], -1)])
import pytest
sp = pytest.importorskip("scipy", minversion="0.7.0")
import numpy as np
import pytensor
from pytensor import sparse
from pytensor.configdefaults import config
from pytensor.sparse.sandbox.sp2 import (
Binomial,
Multinomial,
Poisson,
multinomial,
poisson,
)
from pytensor.tensor.type import lscalar, lvector, scalar
from tests import unittest_tools as utt
from tests.sparse.test_basic import as_sparse_format
class TestPoisson(utt.InferShapeTester):
x = {}
a = {}
for format in sparse.sparse_formats:
variable = getattr(pytensor.sparse, format + "_matrix")
a_val = np.array(
np.random.default_rng(utt.fetch_seed()).integers(1, 4, size=(3, 4)) - 1,
dtype=pytensor.config.floatX,
)
x[format] = variable()
a[format] = as_sparse_format(a_val, format)
def setup_method(self):
super().setup_method()
self.op_class = Poisson
def test_op(self):
for format in sparse.sparse_formats:
f = pytensor.function([self.x[format]], poisson(self.x[format]))
tested = f(self.a[format])
assert tested.format == format
assert tested.dtype == self.a[format].dtype
assert np.allclose(np.floor(tested.data), tested.data)
assert tested.shape == self.a[format].shape
def test_infer_shape(self):
for format in sparse.sparse_formats:
self._compile_and_check(
[self.x[format]],
[poisson(self.x[format])],
[self.a[format]],
self.op_class,
)
class TestBinomial(utt.InferShapeTester):
n = scalar(dtype="int64")
p = scalar()
shape = lvector()
_n = 5
_p = 0.25
_shape = np.asarray([3, 5], dtype="int64")
inputs = [n, p, shape]
_inputs = [_n, _p, _shape]
def setup_method(self):
super().setup_method()
self.op_class = Binomial
def test_op(self):
for sp_format in sparse.sparse_formats:
for o_type in sparse.float_dtypes:
f = pytensor.function(
self.inputs, Binomial(sp_format, o_type)(*self.inputs)
)
tested = f(*self._inputs)
assert tested.shape == tuple(self._shape)
assert tested.format == sp_format
assert tested.dtype == o_type
assert np.allclose(np.floor(tested.todense()), tested.todense())
def test_infer_shape(self):
for sp_format in sparse.sparse_formats:
for o_type in sparse.float_dtypes:
self._compile_and_check(
self.inputs,
[Binomial(sp_format, o_type)(*self.inputs)],
self._inputs,
self.op_class,
)
class TestMultinomial(utt.InferShapeTester):
p = sparse.csr_matrix()
_p = sp.sparse.csr_matrix(
np.asarray(
[
[0.0, 0.5, 0.0, 0.5],
[0.1, 0.2, 0.3, 0.4],
[0.0, 1.0, 0.0, 0.0],
[0.3, 0.3, 0.0, 0.4],
],
dtype=config.floatX,
)
)
def setup_method(self):
super().setup_method()
self.op_class = Multinomial
def test_op(self):
n = lscalar()
f = pytensor.function([self.p, n], multinomial(n, self.p))
_n = 5
tested = f(self._p, _n)
assert tested.shape == self._p.shape
assert np.allclose(np.floor(tested.todense()), tested.todense())
assert tested[2, 1] == _n
n = lvector()
f = pytensor.function([self.p, n], multinomial(n, self.p))
_n = np.asarray([1, 2, 3, 4], dtype="int64")
tested = f(self._p, _n)
assert tested.shape == self._p.shape
assert np.allclose(np.floor(tested.todense()), tested.todense())
assert tested[2, 1] == _n[2]
def test_infer_shape(self):
self._compile_and_check(
[self.p], [multinomial(5, self.p)], [self._p], self.op_class, warn=False
)
Markdown 格式
0%
您添加了 0 到此讨论。请谨慎行事。
请先完成此评论的编辑!
注册 或者 后发表评论