ND convolution with im2col

[jeffdonahue]

This PR extends convolution to N spatial axes, where Caffe's current convolution supports only 2D convolution (with 2 spatial axes: height and width). For 2D convolution, this implementation doesn't compare favorably with the existing one -- I haven't done much benchmarking, but I believe it's 25-75% slower on both CPU and GPU. So before this could be merged, I'd need to restore the existing implementation and use it as the default "engine" for 2D convolutions (but this more destructive version makes it easier to tell what I was thinking from looking at the diff). If anyone has any suggestions on improving the performance or thoughts on why it might be so much slower, I'd love to hear them.

Edit: benchmarking this on alexnet, it's about 33% slower:

@ master:

I0305 21:07:25.042047 22060 caffe.cpp:271] Average Forward pass: 486.327 ms.
I0305 21:07:25.042064 22060 caffe.cpp:273] Average Backward pass: 824.147 ms.
I0305 21:07:25.042079 22060 caffe.cpp:275] Average Forward-Backward: 1310.68 ms.

@ nd-convolution:

I0305 21:02:03.827594 12909 caffe.cpp:271] Average Forward pass: 681.38 ms.
I0305 21:02:03.827608 12909 caffe.cpp:273] Average Backward pass: 1068.98 ms.
I0305 21:02:03.827623 12909 caffe.cpp:275] Average Forward-Backward: 1750.56 ms.

[barhomi]

Is this using cudnn v2 as a backend for the Nd-Conv? if that's the case, I think nvidia's nd-conv (only 3D for now) is not as tuned as their 2D conv, from the release notes:

"As a BETA preview in this release, the convolution forward, convolution
weight and data gradient, and cudnnSetTensor/cudnnScaleTensor routines
now support 3D datasets through the “Nd” interface. Note, however, that
these 3D routines have not yet been tuned for optimal performance. This will
be improved in future releases."

[jeffdonahue]

No, this doesn't touch CuDNN, it only generalizes the im2col convolution implementation (which predates CuDNN).

[barhomi]

ah, sorry! my bad

On Tue, Mar 31, 2015 at 5:34 PM, Jeff Donahue [email protected]
wrote:

No, this doesn't touch CuDNN, it only generalizes the im2col convolution
implementation (which predates CuDNN).

—
Reply to this email directly or view it on GitHub
#2049 (comment).

Youssef Barhomi
PhD candidate
Brown University and University of Pierre and Marie Curie
T: +1 (617) 797 9929 | GMT -5:00 | http://barhomi.github.io

[wkal]

Could you provide a demo to show how to use it? Otherwise, there will be huge high learning curve to test and use your work. Thanks!

[avalada]

@jeffdonahue Does this also support 1D conv?

[jeffdonahue]

@avalada yes, any N >= 0 should theoretically be supported. In practice, 0D convolution -- scalar multiplication -- probably doesn't work, but should and would make a great unit test. I expect 1-10D convolution to work out of the box with this; >10 won't work on GPU -- you'd have to add your case to the switch statements in 718802e.

Also, 1D convolution is supported by the current implementation as well; just set either the width or height to 1. Theoretically, doing 1D convolution using an ND implementation could/should be more efficient than using a 2D implementation with a singleton dim, but with the apparently large overhead in the 2D case, I would be surprised if that's the case here -- you're probably better off sticking with the existing 2D implementation. (But I'd be very interested to know the comparison if you decide to benchmark.)

[ghost]

Hey Jeff,

Is there any chance you could link to an example prototxt making use of this pull request? It would be nice to have that to get started.

[jmerkow]

I don't think there are changes needed in the prototxt to use this PR. Just set your dims using repeated values in the prototxt. i.e.:

name: "4d-net"
input: "data"
input_dim: 1
input_dim: 3 
input_dim: 5
input_dim: 5
input_dim: 5
input_dim: 5

The channel axis defaults to 1 (in this net there are 3 channels). If you want nd kernels just repeat kernel for each dim, instead of using kernel_h,kernel_w.

The notes in caffe.proto describe it pretty well.

[jeffdonahue]

Thanks @jmerkow -- there's a slight correction as input_dim is deprecated in favor of input_shape, as input_dim only works for 4D blobs. It should be:

name: "4d-net"
input: "data"
input_shape { dim: 1 dim: 3 dim: 5 dim: 5 dim: 5 dim: 5 }

Or a full version with DummyData that you should be able to run (didn't test but it should work, possibly needing minor typo fixing):

name: "4d-net"
# -> data: 10 x 3 x 2 x 3 x 4 x 5
# -> label: 10
layer {
  type: "DummyData"
  top: "data"
  top: "label"
  dummy_data_param {
    shape { dim: 10 dim: 3 dim: 2 dim: 3 dim: 4 dim: 5 }
    shape { dim: 10 }
  }
}
# -> conv-out: 10 x 8 x 1 x 1 x 1 x 1
layer {
  type: "Convolution"
  bottom: "data"
  top: "conv-out"
  convolution_param {
    num_output: 8
    # specifies the index of the "channels" axis --
    # may be omitted as 1 is the default
    axis: 1
    kernel_size: 2
    kernel_size: 3
    kernel_size: 4
    kernel_size: 5  # i could just specify "kernel_size: 2" if i wanted 2x2x2x2 filters
  }
}
layer {
  type: "SoftmaxWithLoss"
  bottom: "conv-out"
  bottom: "label"
  top: "loss"
}

[ghost]

@jeffdonahue thanks for the reference. Here is a debugged version of Jeff's prototxt if anyone else is interested (the layers needed names and the SoftmaxWithLoss layer doesn't like >4d blobs):

name: "4d-net"
# -> data: 10 x 3 x 2 x 3 x 4 x 5
# -> label: 10
layer {
  name: "dummy"
  type: "DummyData"
  top: "data"
  top: "label"
  dummy_data_param {
    shape { dim: 10 dim: 3 dim: 2 dim: 3 dim: 4 dim: 5 }
    shape { dim: 10 dim: 1 }
  }
}
# -> conv-out: 10 x 8 x 1 x 1 x 1 x 1
layer {
  name: "conv-out"
  type: "Convolution"
  bottom: "data"
  top: "conv-out"
  convolution_param {
    num_output: 8
    # specifies the index of the "channels" axis --
    # may be omitted as 1 is the default
    axis: 1
    kernel_size: 2
    kernel_size: 3
    kernel_size: 4
    kernel_size: 5  # i could just specify "kernel_size: 2" if i wanted 2x2x2x2 filters
  }
}
layer {
  name: "ip-out"
  type: "InnerProduct"
  bottom: "conv-out"
  top: "ip-out"
  inner_product_param {
    num_output: 1
  }
}
layer {
  name: "loss"
  type: "SoftmaxWithLoss"
  bottom: "ip-out"
  bottom: "label"
  top: "loss"
}

[jmerkow]

@jeffdonahue , If you change Line 150 in filler.hpp to remove legacy calls (i.e. use blob->shape(0)).
You'd probably want to change the sqrt(3) to account for dimensions too.
You can fill your weights with xavier for nd convolution. There are a number of legacy calls in that file that need to be updated.
Im happy to make a PR to update.

https://github.com/jeffdonahue/caffe/blob/718802ef46045c3efe6ef6b5035702d692138010/include/caffe/filler.hpp#L150

[tomdeschamps]

Thanks for sharing the prototxt @Russell91. I'm trying to use this with ND data (N>=3).
Do you have an idea how I can create a lmdb database from a bunch of ND datasets? The example convert_imageset uses caffe::Datum which is limited to height/weight/channel.

[jmerkow]

@tomdeschamps I would try a hdf5 data layer. I believe those can be used to load ND images.

[tomdeschamps]

Thanks @jmerkow. Yes, I'm trying to load it as ND using the hdf5 data layer, but I have an error in LegacyShape() in blob.hpp:141 "Cannot use legacy accessors on Blobs with > 4 axes".
My HDF5Data dimensions is typically (batch size, 1, Z, Y, X) These are batches of grayscale volumes of size XxYxZ. Any idea what I should do? Did anybody actually tried to load ND data in this branch?

[jeffdonahue]

Sorry for the trouble -- there are indeed a lot of places in the code that still use the legacy blob dim accessors (num, channels, height, width). As a quick fix, you might be able to remove that check, and depending on what layers you're using, things might be fine...

The legacy accessors should be removed from most places, definitely in HDF5DataLayer -- PRs welcome. I'm slightly hesitant to say they should be removed everywhere as they can make the code clearer where they're correct -- I introduced the shape call with ND blobs in #1970, and decided against a global find/replace of the legacy calls for that reason. Maybe they should only remain in the still-4D-only vision layers (PoolingLayer, LRNLayer). Any thoughts from other Caffe devs?

[tomdeschamps]

Yes that seems to work with the check removed. I used @jmerkow nd-pooling branch. However HDF5Data seems to be handled very differently (no scaling, not sure we write whole batches in each .h5 file, etc...). Is there a documentation on how caffe deals with this format?

[dzhwinter]

@jeffdonahue @tomdeschamps tested, check remove can be compile successfully. I'm wondering that nd convolution and nd pooling in @jmerkow can be used together now?

[dzhwinter]

Hi, is it possible to implement an cudnn 3d convolution? Flatten 3d volumn into 2d dim or something like that.
Because it really hurt the speed without support cudnn. On another side, 3d convolution always process videos, which need more compute resources.

[Tgaaly]

@jeffdonahue @Russell91 This is great work! Thanks. Today is my first day brewing Caffe!

I'm trying to use this nd-conv with a "Flatten" layer. The "Flatten" layer does not work on blobs with > 4 axes so it fails when the n-dimensions are larger than 2. Is there a fix/solution for this? Would appreciate any advise/help.

I0625 16:20:14.583077  4782 net.cpp:380]  flatdata <- data_dummy_0_split_1
I0625 16:20:14.583091  4782 net.cpp:338] flatdata -> flatdata
I0625 16:20:14.583103  4782 net.cpp:113] Setting up flatdata
F0625 16:20:14.583122  4782 blob.hpp:141] Check failed: num_axes() <= 4 (5 vs. 4) Cannot use legacy accessors on Blobs with > 4 axes.

[Tgaaly]

I disabled the check on line 141 of blob.hpp (shown below) and it ran. Will that cause any problems? Actually looking back at the thread - this seems to be the consensus of others and it make it work.

CHECK_LE(num_axes(), 4)<< "Cannot use legacy accessors on Blobs with > 4 axes.";

[ToruHironaka]

@tomdeschamps, @jmerkow, you guys mentioned about loading ND images with hdf5 earlier. I try to load 3D volume images (Width, height, and Depth) too. I want to make sure a couple of things. I followed your conversations above. I do not think I can convert my 3D image files into lmdb or leveldb. Am I right? It looks like hdf5 data layer is only the way to load my 3D images to ND blob. Have you successfully loaded 3D images? If so, please give me some advise how to load 3D images.

[tomdeschamps]

I think the nd-pooling branch is based on the nd-convolution branch. See
its description.

On Sat, Jun 13, 2015 at 5:51 AM, dzhwinter [email protected] wrote:

@jeffdonahue https://github.com/jeffdonahue @tomdeschamps
https://github.com/tomdeschamps tested, check remove can be compile
successfully. I'm wondering that nd convolution and nd pooling in @jmerkow
https://github.com/jmerkow can be used together now?

—
Reply to this email directly or view it on GitHub
#2049 (comment).

[ToruHironaka]

@Tgaaly, I assume that you stacked a series of 2D images to make 3D image dataset in hdf5 and list hdf5 files in train.txt and test.txt for each classes. Is it so? Sorry, I have been working on something else. It's been for awhile for getting in touch with you.

[ToruHironaka]

Hi, All

I converted 3D images (2 CT-scanned human brain images) into HDF5 datasets with (Number of images, Width, Length, channel [not specified for grayscale]) in python. I created HDF5 dataset files for each 3D volume images and listed these HDF5 files on train.txt and test.txt. Then, I defined a net with below code from caffe/examples/02-brewing-logreg.ipynb

Note: below codes are the original codes from the example so I modified this code for my dataset and net.

def logreg(hdf5, batch_size):
# logistic regression: data, matrix multiplication, and 2-class softmax loss
n = caffe.NetSpec()
n.data, n.label = L.HDF5Data(batch_size=batch_size, source=hdf5, ntop=2)
n.ip1 = L.InnerProduct(n.data, num_output=2, weight_filler=dict(type='xavier'))
n.accuracy = L.Accuracy(n.ip1, n.label)
n.loss = L.SoftmaxWithLoss(n.ip1, n.label)
return n.to_proto()

with open('examples/hdf5_classification/logreg_auto_train.prototxt', 'w') as f:
f.write(str(logreg('examples/hdf5_classification/data/train.txt', 10)))

with open('examples/hdf5_classification/logreg_auto_test.prototxt', 'w') as f:
f.write(str(logreg('examples/hdf5_classification/data/test.txt', 10)))

After that, I run test with below code from the same example

caffe.set_mode_gpu()
solver = caffe.get_solver('examples/hdf5_classification/solver.prototxt')
solver.solve()

accuracy = 0
batch_size = solver.test_nets[0].blobs['data'].num
test_iters = int(len(Xt) / batch_size)
for i in range(test_iters):
solver.test_nets[0].forward()
accuracy += solver.test_nets[0].blobs['accuracy'].data
accuracy /= test_iters

print("Accuracy: {:.3f}".format(accuracy))

I got the results like below and it looks alright. Can anyone confirm my way of building caffe 3D (Depth,Width,Height,Channel [channel ignored if gray scale image]) model is correct or not?

Results:

I1120 12:43:45.839637 28983 solver.cpp:734] Snapshotting solver state to binary proto file/hdf_FT_iter_10000.solverstate
I1120 12:43:45.846904 28983 solver.cpp:326] Iteration 10000, loss = 0
I1120 12:43:45.847026 28983 solver.cpp:346] Iteration 10000, Testing net (#0)
I1120 12:43:45.847059 28983 net.cpp:781] Copying source layer data
I1120 12:43:45.847087 28983 net.cpp:781] Copying source layer label_data_1_split
I1120 12:43:45.847115 28983 net.cpp:781] Copying source layer ip1
I1120 12:43:45.847146 28983 net.cpp:781] Copying source layer ip1_ip1_0_split
I1120 12:43:45.847172 28983 net.cpp:781] Copying source layer accuracy
I1120 12:43:45.847198 28983 net.cpp:781] Copying source layer loss
I1120 12:43:45.847236 28983 hdf5_data_layer.cu:33] Looping around to first file.
I1120 12:43:45.847262 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTBraintest.h5
I1120 12:43:45.856613 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:45.862017 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTImagetest.h5
I1120 12:43:45.874867 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:45.879843 28983 hdf5_data_layer.cu:33] Looping around to first file.
I1120 12:43:45.879860 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTBraintest.h5
.
.
I1120 12:43:46.574106 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:46.579458 28983 solver.cpp:414] Test net output #0: accuracy = 0.9
I1120 12:43:46.579494 28983 solver.cpp:414] Test net output #1: loss = 8.73365 (* 1 = 8.73365 loss)
I1120 12:43:46.579507 28983 solver.cpp:331] Optimization Done.
I1120 12:43:46.580368 28983 hdf5_data_layer.cu:33] Looping around to first file.
I1120 12:43:46.580379 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTBraintest.h5
I1120 12:43:46.589470 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:46.598145 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTImagetest.h5
I1120 12:43:46.606251 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:46.615324 28983 hdf5_data_layer.cu:33] Looping around to first file.
I1120 12:43:46.615339 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTBraintest.h5
I1120 12:43:46.624083 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:46.632520 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTImagetest.h5
I1120 12:43:46.641490 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
I1120 12:43:46.649895 28983 hdf5_data_layer.cu:33] Looping around to first file.
I1120 12:43:46.649909 28983 hdf5_data_layer.cpp:29] Loading HDF5 file: /CTBraintest.h5
I1120 12:43:46.657821 28983 hdf5_data_layer.cpp:68] Successully loaded 20 rows
Accuracy: 0.880

[jeffdonahue]

(Number of images, Width, Length, channel [not specified for grayscale])

The 'channels' axis should be right after the batch axis, so the shape should be (batch_size, channels, *spatial_axes). The ordering of spatial axes is arbitrary. (The only reason channels are treated separately from spatial axes is to support groups, which make the layer act like a locally-connected layer along the channel axis.)

[ToruHironaka]

@jeffdonahue, Thanks for answering me my question, You stated earlier in this pull. Sorry, I missed it. I set my dimensions as follow: depth (stuck of 3D volume), channel (RGB), width, length. Caffe run okay and the result was better. I was probably training width and channel instead of height in my last run. Just one more question, can I use the same method for multi-channel or spectrum images in hdf5?

[jmerkow]

You can load whatever you want with hdf5 as long as its sticks to the NxCxSxSx....(where S is a spatial dim). You can have multiple channels in each image, or multiple batches in each file. I typically stick to a batch of 1, and increase/decrease with the batch_size param. But I don't think you need to, for example if you want images grouped into pre-determined batches.
You can also put multiple images into each hdf5 layer and reference them with different tops. These can be concatenated/silenced if you want to try different combinations.
--Jameson

[ToruHironaka]

@jeffdonahue @jmerkow @Tgaaly, I previously forgot to set convolution_param so my last train was not 3D. It was 2D because I did not add kernel_size. I referred Tgaaly's model layer above but I got the error below.

[libprotobuf ERROR google/protobuf/text_format.cc:245] Error parsing text-format caffe.NetParameter: 32:16: Non-repeated field "kernel_size" is specified multiple times.

I think my hdf5 dataset's dimension is wrong. I have my dataset is (number of files, channel, width, length). I think I am suppose set my dataset dinmesion (batch_size, channel, width, length, depth). Then, add 3 kernel_size in convolution_param.

Do you have any suggestion?

[naibaf7]

@ToruHironaka
If you are interested, I am currently working on a python interface for custom network training/processing in 3D.
It should simplify dataset loading because you can manipulate/preprocess in python and then load it to Caffe via MemoryDataLayers.

Branch: https://github.com/naibaf7/caffe (includes OpenCL support and 3D max pooling support)
Interface: https://github.com/naibaf7/PyGreentea

it's still work in progress but if you can provide me your network and how the data is formatted I might be able to prepare a working python script for you.

[ToruHironaka]

@naibaf7
Thanks, I will check your Branch out.

[anewlearner]

@ToruHironaka
I’m new to caffe and have just run the examples of mnist and cifar10. Now I want to train my own network using 3D CT slices base on caffe-windows. The data I have is in the form of mat(N_C_Sx_Sy_Sz, S:spatial). I read the discussion above, but still confused.
The process to train 2D image is as follows.

1. Changing the data format into Leveldb or LMDB(convert_imageset.cpp).
2. Mean(compute_image_mean.cpp).
3. Train(train_net.cpp).
   My understanding of train 3D image is as follows. I’m not sure whether is correct.
4. Using jeffdonahue’s project(changing caffe.proto,convolution and blobs).
5. Changing the data format into hdf5(I have the original format of mat,so no change).
6. Mean(the input of compute_image_mean should be leveldb,so do I need to write a new function to count the mean?).
7. Train(Can I use train_net.cpp here).

Does anybody has a examples of this?
Thanks!

Ps: I don’t know anything about python.

[ToruHironaka]

@anewlearner

1. Using jeffdonahue’s project(changing caffe.proto,convolution and blobs). Yes
2. Changing the data format into hdf5(I have the original format of mat,so no change). Yes
3. Mean(the input of compute_image_mean should be leveldb,so do I need to write a new function to count the mean?). I did not use mean.
4. Train(Can I use train_net.cpp here). I could use caffe train command to train with this promotion of caffe

You have to use hdf data format for 3D-convolution in this promotion of caffe. I wrote a python script to convert my CT image files (Width, Height, Depth) into hdf data file. Then, I could train my 3D hdf datasets with this promotion of caffe. It worked but I did not get good results yet. My accuracy was very low like 0.4~0.6 and loss was always high like 1.5 or 1.6.

I am now troubleshooting my image-to-hdf python scripts. I tested my python script to create 2D dataset and trained in the official caffe. I got the accuracy about 0.87 and loss was about 0.62. Then, I used other person's image-to-hdf matlab script to create hdf datasets with the same images and trained exactly the same way as my python script test. It got accuracy about 0.88 and loss was about 0.2. I created lmdb datasets with the same images by using caffe conversion command, which you used for 2D images. I got the accuracty about 0.93 and loss was 0.35. So, my image-to-hdf python conversion script was obviously worst. I am finding my data conversion problem now.

This promotion of caffe accepted 3D dataset in hdf5 and it worked. Also, many people confirmed it. You should try it out. If you need my help, let me know because it helps my data conversion problem. If your hdf datasets work, you got my answer.

[anewlearner]

@ToruHironaka
Thanks for your answer!
I plan to change the data with the matlab script store2hdf5 provided by caffe.
The process of mean may help the net to converge faster, you can try it later.

[SiddhantRanade]

Has anyone successfully gotten ND-Pooling to work? ND Convolution works without issues (from the master branch of BVLC-Caffe)

[naibaf7]

@pietromaximoff
https://github.com/BVLC/caffe/tree/opencl
in this branch I have implemented MaxPooling for ND, and it works successfully.
Would be great though if someone has the time to complete this for the other pooling types & make a PR to master.
Feel free to use the code there as a starting point.

Code is here:
https://github.com/BVLC/caffe/blob/opencl/src/caffe/layers/pooling_layer.cu
(MaxPoolingNDForward/MaxPoolingNDBackward)
and the Reshape + constructor from here:
https://github.com/BVLC/caffe/blob/opencl/src/caffe/layers/pooling_layer.cpp

[aliimran9010]

@jeffdonahue Hi, i am new to caffe. I used the "nd-convolution" branch. It gives me the
Check failed: num_kernel_dims == 1 || num_kernel_dims == num_spatial_axes_ kernel_size must be specified once, or once per spatial dimension (kernel_size specified 3 times; 2 spatial dims).

How can i resolve it

[SiddhantRanade]

@naibaf7 I get the following error when I try to use the opencl branch of caffe (with python):

pycaffe.py:13: RuntimeWarning: to-Python converter for std::vector<int, std::allocator > already registered; second conversion method ignored.
from ._caffe import
Floating point exception

I have no idea whatsoever what this error is and how to resolve it. Do you know how I can fix this?

[paulcx]

@ToruHironaka Do you have an example of how to train such data format of CT images on Caffe with 3D convolution?

[ToruHironaka]

@paulcx I wrote a python script for converting image files into hdf5 format and I followed this promotion's thread above. I could trained models but I did not get good results so I did something wrong.

[xjtuljy]

@ToruHironaka Is hdf5 the only format that work with this PR? Did you try N-D max pooling together with N-D convolution?

[ToruHironaka]

@xjtuljy yes, hdf5 is the only format for this promotion. I tried to train my 3D-CNN for ND-Pooling with Promotion 2442 and 2824. They ran and seemed to be working but my result was bad, so I think I am doing something wrong with my training.