A simple CNN for the MNIST datasets – I – CNN basics

Posted on 25. May 2020 by eremo

In a previous article series
A simple Python program for an ANN to cover the MNIST dataset – I – a starting point
we have played with a Python/Numpy code, which created a configurable and trainable “Multilayer Perceptron” [MLP] for us. See also
MLP, Numpy, TF2 – performance issues – Step I – float32, reduction of back propagation
for ongoing code and performance optimization.

A MLP program is useful to study multiple topics in Machine Learning [ML] on a basic level. However, MLPs with dense layers are certainly not at the forefront of ML technology – though they still are fundamental bricks in other more complicated architectures of “Artifical Neural Networks” [ANNs]. During my MLP experiments I became sufficiently acquainted with Python, Jupyter and matplotlib to make some curious first steps into another field of Machine Learning [ML] now: “Convolutional Neural Networks” [CNNs].

CNNs on my level as an interested IT-affine person are most of all fun. Nevertheless, I quickly found out that a somewhat systematic approach is helpful – especially if you later on want to use the Tensorflow’s API and not only Keras. When I now write about some experiments I did and do I summarize my own biased insights and sometimes surprises. Probably there are other hobbyists as me out there who also fight with elementary points in the literature and practical experiments. Books alone are not enough … I hope to deliver some practical hints for this audience. The present articles are, however, NOT intended for ML and CNN experts. Experts will almost certainly not find anything new here.

Although I address CNN-beginners I assume that people who stumble across this article and want to follow me through some experiments have read something about CNNs already. You should know fundamentals about filters, strides and the basic principles of convolution. I shall comment on all these points but I shall not repeat the very basics. I recommend to read relevant chapters in one of the books I recommend at the end of this article. You should in addition have some knowledge regarding the basic structure and functionality of a MLP as well as “gradient descent” as an optimization technique.

The objective of this introductory mini-series is to build a first simple CNN, to apply it to the MNIST dataset and to visualize some of the elementary “features” a CNN allegedly detects in the images of handwritten digits – at least according to many authors in the field of AI. We shall use Keras (with the Tensorflow 2.2 backend and CUDA 10.2) for this purpose. And, of course, a bit of matplotlib and Python/Numpy, too. We are working with MNIST images in the first place – although CNNs can be used to analyze other types of input data. After we have covered the simple standard MNIST image set, we shall also work a bit with the so called “MNIST fashion” set.

But in this article I start with some introductory words on the structure of CNNs and the task of its layers. We shall use the information later on as a reference. In the second article we shall set up and test a simple version of a CNN. Further articles will then concentrate on visualizing what a trained CNN reacts to and how it modifies and analyzes the input data on its layers.

Why CNNs?

When we studied an MLP in combination with the basic MNIST dataset of handwritten digits we found that we got an improvement in accuracy (for the same setup of dense layers) when we pre-processed the data to find “clusters” in the image data before training. Such a process corresponds to detecting parts of an MNIST image with certain gray-white pixel constellations. We used Scikit-Learn’s “MiniBatchKMeans” for this purpose.

We saw that the identification of 40 to 70 cluster areas in the images helped the MLP algorithm to analyze the MNIST data faster and better than before. Obviously, training the MLP with respect to combinations of characteristic sub-structures of the different images helped us to classify them as representations of digits. This leads directly to the following question:

What if we could combine the detection of sub-structures in an image with the training process of an ANN?

CNNs seem to the answer! According to teaching books they have the following abilities: They are designed to detect elementary structures or patterns in image data (and other data) systematically. In addition they are enabled to learn something about characteristic compositions of such elementary features during training. I.e., they detect more abstract and composite features specific for the appearance of certain objects within an image. We speak of a “feature hierarchy“, which a CNN can somehow grasp and use – e.g. for classification tasks.

While a MLP must learn about pixel constellations and their relations on the whole image area, CNNs are much more flexible and even reusable. They identify and remember elementary sub-structures independent of the exact position of such features within an image. They furthermore learn “abstract concepts” about depicted objects via identifying characteristic and complex composite features on a higher level.

This simplified description of the astonishing capabilities of a CNN indicates that its training and learning is basically a two-fold process:

Detecting elementary structures in an image (or other structured data sets) by filtering and extracting patterns within relatively small image areas. We shall call these areas “filter areas”.
Constructing abstract characteristic features out of the elementary filtered structural elements. This corresponds to building a “hierarchy” of significant features for the classification of images or of distinguished objects or of the positions of such objects within an image.

Now, if you think about the MNIST digit data we understand intuitively that written digits represent some abstract concepts like certain combinations of straight vertical and horizontal line elements, bows and line crossings. The recognition of certain feature combinations of such elementary structures would of course be helpful to recognize and classify written digits better – especially when the recognition of the combination of such features is independent of their exact position on an image.

So, CNNs seem to open up a world of wonders! Some authors of books on CNNs, GANs etc. praise the ability to react to “features” by describing them as humanly interpretable entities as e.g. “eyes”, “feathers”, “lips”, “line segments”, etc. – i.e. in the sense of entity conceptions. Well, we shall critically review this idea, which I think is a misleading over-interpretation of the capacities of CNNs.

Filters, kernels and feature maps

An important concept behind CNNs is the systematic application of (various) filters (described and defined by so called “kernels”).

A “filter” defines a kind of masking pixel area of limited small size (e.g. 3×3 pixels). A filter combines weighted output values at neighboring nodes of a input layer in a specific defined way. It processes the offered information in a defined area always in the same fixed way –
independent of where the filter area is exactly placed on the (bigger) image (or a processed version of it). We call a processed version of an image a “map“.

A specific type of CNN layer, called a “Convolution Layer” [Conv layer], and a related operational algorithm let a series of such small masking areas cover the complete surface of an image (or a map). The first Conv layer of a CNN filters the information of the original image information via a multitude of such masking areas. The masks can be arranged overlapping, i.e. they can be shifted against each other by some distance along their axes. Think of the masking filter areas as a bunch of overlapping tiles covering the image. The shift is called stride.

The “filter” mechanism (better: the mathematical recipe) of a specific filter remains the same for all of its small masking areas covering the image. A specific filter emphasizes certain parts of the original information and suppresses other parts in a defined way. If you combine the information of all masks you get a new (filtered) representation of the image – we speak of a “feature map” – sometimes with a smaller size than the original image (or map) the filter is applied to. The blending of the original data with a filtering mask creates a “feature map“, i.e. a filtered view onto the input data. The blending process is called “convolution” (due to the related mathematical operations).

The picture below sketches the basic principle of a 3×3-filter which is applied with a constant stride of 2 along each axis of the image:

Convolution is not so complicated as it sounds. It means: You multiply the original data values in the covered small area by factors defined in the filter’s kernel and add the resulting values up to get a a distinct value at a defined position inside the map. In the given example with a stride of 2 we get a resulting feature map of 4×4 out of a original 9×9 (image or map).

Note that a filter need not be defined as a square. It can have a rectangular (n x m) shape with (n, m) being integers. (In principle we could also think of other tile forms as e.g. hexagons – as long as they can seamlessly cover a defined plane. Interesting, although I have not seen a hexagon based CNN in the literature, yet).

A filter’s kernel defines factors used in the convolution operation – one for each of the (n x m) defined points in the filter area.

Note also that filters may have a “depth” property when they shall be applied to three-dimensional data sets; we may need a depth when we cover colored images (which require 3 input layers). But let us keep to flat filters in this introductory discussion …

Now we come to a central question: Does a CNN Conv layer use just one filter? The answer is: No!

A Conv layer of a CNN you allows for the construction of multiple different filters. Thus we have to deal with a whole bunch of filters per each convolutional layer. E.g. 32 filters for the first convolutional layer and 64 for the second and 128 for the third. The outcome of the respective filter operations is the creation is of equally many “feature maps” (one for each filter) per convolutional layer. With 32 different filters on a Conv layer we would thus build 32 maps at this layer.

This means: A Conv layer has a multitude of sub-layers, i.e. “maps” which result of the application of different filters on previous image or map data.

You may have guessed already that the next step of abstraction is:
You can apply filters also to “maps” of previous filters, i.e. you can chain convolutions. Thus, feature maps are either connected to the image (1st Conv layer) or to the maps of a previous layer.

By using a sequence of multiple Conv layers you cover growing areas of the original image. Everything clear? Probably not …

Filters and their related weights are the end products of the training and optimization of a CNN!

When I first was confronted with the concept of filters, I got confused because many authors only describe the basic technical details of the “convolution” mechanism. They explain with many words how a filter and its kernel work when the filtering area is “moved” across the surface of an image. They give you pretty concrete filter examples; very popular are straight lines and crosses indicated by “ones” as factors in the filter’s kernel and zeros otherwise. And then you get an additional lecture on strides and padding. You have certainly read various related passages in books about ML and/or CNNs. A pretty good example for this “explanation” is the (otherwise interesting and helpful!) book of Deru and Ndiaye (see the bottom of this article. I refer to the introductory chapter 3.5.1 on CNN architectures.)

Well, the technical procedure is pretty easy to understand from drawings as given above – the real question that nags in your brain is:

“Where the hell do all the different filter definitions come from?”

What many authors forget is a central introductory sentence for beginners:

A filter is not given a priori. Filters (and their kernels) are systematically constructed and build up during the training of a CNN; filters are the end products of a learning and optimization process every CNN must absolve.

This means: For a given problem or dataset you do not know in advance what the “filters” (and their defining kernels) will look like after training (aside of their pixel dimensions already fixed by the CNN’s layer definitions). The “factors” of a filter used in the convolution operation are actually weights, whose final values are the outcome of a learning process. Just as in MLPs …

Nothing is really “moved” …

Another critical point is the somewhat misleading analogy of “moving” a filter across an image’s or map’s pixel surface. Nothing is ever actually “moved” in a CNN’s algorithm. All masks are already in place when the convolution operations are performed:

Every element of a specific e.g. 3×3 kernel corresponds to “factors” for the convolution operation. What are these factors? Again: They are nothing else but weights – in exactly the same sense as we used them in MLPs. A filter kernel represents a set of weight-values to be multiplied with original output values at the “nodes” in other layers or maps feeding input to the nodes of the present map.

Things become much clearer if you imagine a feature map as a bunch of arranged “nodes”. Each node of a map is connected to (n x m) nodes of a previous set of nodes on a map or layer delivering input to the Conv layer’s maps.

Let us look at an example. The following drawing shows the connections from “nodes” of a feature map “m” of a Conv layer L_(N+1) to nodes of two different maps “1” and “2” of
Conv layer L_N. The stride for the kernels is assumed to be just 1.

In the example the related weights are described by two different (3×3) kernels. Note, however, that each node of a specific map uses the same weights for connections to another specific map or sub-layer of the previous (input) layer. This explains the total number of weights between two sequential Conv layers – one with 32 maps and the next with 64 maps – as (64 x 32 x 9) + 64 = 18496. The 64 extra weights account for bias values per map on layer L_(N+1). (As all nodes of a map use fixed bunches of weights, we only need exactly one bias value per map).

Note also that a stride is defined for the whole layer and not per map. Thus we enforce the same size of all maps in a layer. The convolutions between a distinct map and all maps of the previous layer L_N can be thought as operations performed on a column of stacked filter areas at the same position – one above the other across all maps of L_N. See the illustration below:

The weights of a specific kernel work together as an ensemble: They condense the original 3×3 pixel information in the filtered area of the connected input layer or a map to a value at one node of the filter specific feature map. Please note that there is a bias weight in addition for every map; however, at all masking areas of a specific filter the very same 9 weights are applied. See the next drawing for an illustration of the weight application in our example for fictitious node and kernel values.

A CNN learns the appropriate weights (= the filter definitions) for a given bunch of images via training and is guided by the optimization of a loss function. You know these concepts already from MLPs …

The difference is that the ANN now learns about appropriate “weight ensembles” – eventually (!) working together as a defined convolutional filter between different maps of neighboring Conv (and/or sampling ) Layers. (For sampling see a separate paragraph below.)

The next picture illustrates the column like convolution of information across the identically positioned filter areas across multiple maps of a previous convolution layer:

The fact that the weight ensemble of a specific filter between maps is always the same, explains, by the way, the relatively (!) small number of weight parameters in deep CNNS.

Intermediate summary: The weights, which represent the factors used by a specific filter operation called convolution, are defined during a training process. The filter, its kernel and the respective weight values are the outcome of a mathematical optimization process – mostly guided by gradient descent.

Activation functions

As in MLPs each Conv layer has an associated “activation function” which is applied at each node of all maps after the resulting values of the convolution have been calculated as
the nodes input. The output then feeds the connections to the next layer. In CNNs for image handling often “Relu” or “Selu” are used as activation functions – and not “sigmoid” which we applied in the MLP code discussed in another article series of this blog.

Tensors

The above drawings indicate already that we need to arrange the data (of an image) and also the resulting map data in an organized way to be able to apply the required convolutional multiplications and summations the right way.

An colored image is basically a regular 3 dimensional structure with a width “w” (number of pixels along the x-axis), a height “h” (number of pixels along the y-axis) and a (color) depth “d” (d=3 for RGB colors).

If you represent the color value at each pixel and RGB-layer by a float you get a bunch of w x h x d float values which we can organize and index in a 3 dimensional Numpy array. Mathematically such well organized arrays with a defined number of axes (rank), a set of numbers describing the dimension along each axis (shape), a data-type, possible operations (and invariance aspects) define an abstract object called a “tensor“. Colored image data can be arranged in 3-dimensional tensors; gray colored images in a pseudo 3D-tensor which has a shape of (n, m, 1). (Keras and Tensorflow want to get imagedata in form of 2D tensors).

Now the important point is: The output data of Conv-layers and their feature maps also represent tensors. A bunch of 32 maps with a defined width and height defines data of a 3D-tensor.

You can imagine each value of such a tensor as the input or output given at a specific node in a layer with a 3-dimensional sub-structure. (In other even more complex data structures than images we would other multi-dimensional data structures.) The weights of a filter kernel describe the connections of the nodes of a feature map on a layer L_N to a specific map of a previous layer. Weights, actually, also define elements of a tensor.

The forward- and backward-propagation operations performed throughout such a complex net during training thus correspond to sequences of tensor-operations – i.e. generalized versions of the np.dot()-product we got to know in MLPs.

You understood already that e.g strides are important. But you do not need to care about details – Keras and Tensorflow will do the job for you! If you want to read a bit look a the documentation of the TF function “tf.nn.conv2d()”.

When we later on train with mini-batches of input data (i.e. batches of images) we get yet another dimension of our tensors. This batch dimension can – quite similar to MLPs – be used to optimize the tensor operations in a vectorized way. See my series on MLPs.

Chained convolutions cover growing areas of the original image

Two sections above I characterized the training of a CNN as a two-fold procedure. From the first drawing it is relatively easy to understand how we get to grasp tiny sub-structures of an image: Just use filters with small kernel sizes!

Fine, but there is probably a second question already arising in your mind:

By what mechanism does a CNN find or recognize a hierarchy of features?

One part of the answer is: Chain convolutions!

Let us assume a first convolutional layer with filters having a stride of 1 and a (3×3) kernel. We get maps with a shape of (26, 26) on this layer. The next Conv layer shall use a (4×4) kernel and also a stride of 1; then we get maps with a shape of (23, 23). A node on the second layer covers (6×6)-arrays on the original image. Two neighboring nodes a total area of (7×7). The individual (6×6)-areas of course overlap.

With a stride of 2 on each Conv-layer the corresponding areas on the original image are (7×7) and (11×11).

So a stack of consecutive (sequential) Conv-layers covers growing areas on the original image. This supports the detection of a pattern or feature hierarchy in the data of the input images.

However: Small strides require a relatively big number of sequential Conv-layers (for 3×3 kernels and stride 2) at least 13 layers to eventually cover the full image area.

Even if we would not enlarge the number of maps beyond 128 with growing layer number, we would get

(32 x 9 + 32) + (64 x 32 +64) + (128 x 64 + 128) + 10 x (128 x 128 + 128) = 320 + 18496 + 73856 + 10*147584 = 1.568 million weight parameters

to take care of!

This number has to be multiplied by the number of images in a mini-batch – e.g. 500. And – as we know from MLPs we have to keep all intermediate output results in RAM to accelerate the BW propagation for the determination of gradients. Too many data and parameters for the analysis of small 28×28 images!

Big strides, however, would affect the spatial resolution of the first layers in a CNN. What is the way out?

Sub-sampling is necessary!

The famous VGG16 CNN uses pairs and triples of convolution chains in its architecture. How does such a network get control over the number of weight parameters and the RAM requirement for all the output data at all the layers?

To get information in the sense of a feature hierarchy the CNN clearly should not look at details and related small sub-fields of the image, only. It must cover step-wise growing (!) areas of the original image, too. How do we combine these seemingly contradictory objectives in one training algorithm which does not lead to an exploding number of parameters, RAM and CPU time? Well, guys, this is the point where we should pay due respect to all the creative inventors of CNNs:

The answer is: We must accumulate or sample information across larger image or map areas. This is the (underestimated?) task of pooling– or sampling-layers.

For me it was just another confusing point in the beginning – until one grasps the real magic behind it. At first sight a layer like a typical “maxpooling” layer seems to reduce information, only; see the next picture:

The drawing explains that we “sample” the information over multiple pixels e.g. by

either calculating an average over pixels (or map node values)
or by just picking the maximum value of pixels or map node values (thereby stressing the most important information)

in a certain defined sub-area of an image or map.

The shift or stride used as a default in a pooling layer is exactly the side length of the pooling area. We thus cover the image by adjacent, non-overlapping tiles! This leads to a substantial decrease of the dimensions of the resulting map! With a (2×2) pooling size by a an effective factor of 2. (You can change the default pooling stride – but think about the consequences!)

Of course, averaging or picking a max value corresponds to information reduction.

However: What the CNN really also will do in a subsequent Conv layer is to invest in further weights for the combination of information (features) in and of substantially larger areas of the original image! Pooling followed by an additional convolution obviously supports hierarchy building of information on different scales of image areas!

After we first have concentrated on small scale features (like with a magnifying glass) we now – in a figurative sense – make a step backwards and look at larger scales of the image again.

The trick is to evaluate large scale information by sampling layers in addition to the small scale information information already extracted by the previous convolutions. Yes, we drop resolution information – but by introducing a suitable mix of convolutions and sampling layers we also force the network systematically to concentrate on combined large scale features, which in the end are really important for the image classification as a whole!

As sampling counterbalances an explosion of parameters we can invest into a growing number of feature maps with growing scales of covered image areas. I.e. we add more and new filters reacting to combinations of larger scale information.

Look at the second to last illustration: Assume that the 32 maps on layer L_N depicted there are the result of a sampling operation. The next convolution gathers new knowledge about more, namely 64 different combinations of filtered structures over a whole vertical stack of small filter areas located at the same position on the 32 maps of layer N. The new information is in the course of training conserved into 64 weight ensembles for 64 maps on layer N+1.

Resulting options for architectures

We can think of multiple ways of combining Conv layers and pooling layers. A simple recipe for small images could be

Layer 0: Input layer (tensor of original image data, 3 color layers or one gray layer)
Layer 1: Conv layer (small 3×3 kernel, stride 1, 32 filters, 32 maps (26×26), analyzes 3×3 overlapping areas)
Layer 2: Pooling layer (2×2 max pooling => 32 (13×13) maps,
a node covers 4×4 non overlapping areas per node on the original image)
Layer 3: Conv layer (3×3 kernel, stride 1, 64 filters, 64 maps (11×11),
a node covers 8×8 overlapping areas on the original image (total effective stride 2))
Layer 4: Pooling layer (2×2 max pooling => 64 maps (5×5),
a node covers 10×10 areas per node on the original image (total effective stride 5), some border info lost)
Layer 5: Conv layer (3×3 kernel, stride 1, 64 filters, 64 maps (3×3),
a node covers 18×18 per node (effective stride 5), some border info lost )

The following picture illustrates the resulting successive combinations of nodes along one axis of a 28×28 image.

Note that I only indicated the connections to border nodes of the Conv filter areas.

The kernel size decides on the smallest structures we look at – especially via the first convolution. The sampling decides on the sequence of steadily growing areas which we then analyze for specific combinations of smaller structures.

Again: It is most of all the (down-) sampling which allows for an effective hierarchical information building over growing larger image areas! Actually we do not really drop information by sampling – instead we give the network a chance to collect and code new information on a higher, more abstract level (via a whole bunch of numerous new weights).

The big advantages of the sampling layers get obvious:

They reduce the numbers of required weights
They reduce the amount of required memory – not only for weights but also for the output data, which must be saved for every layer, map and node.
They reduce the CPU load for FW and BW propagation
They also limit the risk of overfitting as some detail information is dropped.

Of course there are many other sequences of layers one could think about. E.g., we could combine 2 to 3 Conv layers before we apply a pooling layer. Such a layer sequence is characteristic of the VGG nets.

Further aspects

Just as MLPs a CNN represents an acyclic graph, where the maps contain increasingly fewer nodes but where the number of maps per layer increases on average.

Questions and objectives for this article series

An interesting question, which seldom is answered in introductory books, is whether two totally independent training runs for a given CNN-architecture applied on the same input data will produce the same filters in the same order. We shall investigate this point in the forthcoming articles.

Another interesting point is: What does a CNN see at which convolution layer?
And even more important: What do the “features” (= basic structural elements) in an image which trigger/activate a specific filter or map, look like?

If we could look into the output at some maps we could possibly see what filters do with the original image. And if we found a way to construct a structured image which triggers a specific filter then we could better understand what patterns the CNN reacts to. Are these patterns really “features” in the sense of conceptual entities? Examples for these different types of visualizations with respect to convolution in a CNN are objectives of this article series.

Conclusion

Today we covered a lot of “theory” on some aspects of CNNs. But we have a sufficiently solid basis regarding the structure and architecture now.

CNNs obviously have a much more complex structure than MLPs: They are deep in the sense of many sequential layers. And each convolutional layer has a complex structure in form of many parallel sub-layers (feature maps) itself. Feature maps are associated with filters, whose parameters (weights) get learned during the training. A map results from covering the original image or a map of a previous layer with small (overlapping) tiles of small filtering areas.

A mix of convolution and pooling layers allows for a look at detail patterns of the image in small areas in lower layers, whilst later layers can focus on feature combinations of larger image areas. The involved filters thus allow for the “awareness” of a hierarchy of features with translational invariance.

Pooling layers are important because they help to control the amount of weight parameters – and they enhance the effectiveness of detecting the most important feature correlations on larger image scales.

All nice and convincing – but the attentive reader will ask: Where and how do we do the classification?
Try to answer this question yourself first.

In the next article we shall build a concrete CNN and apply it to the MNIST dataset of images of handwritten digits. And whilst we do it I deliver the answer to the question posed above. Stay tuned …

Literature

“Advanced Machine Learning with Python”, John Hearty, 2016, Packt Publishing – See chapter 4.

“Deep Learning mit Python und Keras”, Francois Chollet, 2018, mitp Verlag – See chapter 5.

“Hands-On Machine learning with SciKit-Learn, Keras & Tensorflow”, 2nd edition, Aurelien Geron, 2019, O’Reilly – See chapter 14.

“Deep Learning mit Tensorflow, keras und Tensorflow.js”, Matthieu Deru, Alassane Ndiaye, 2019, Rheinwerk Verlag, Bonn – see chapter 3

Further articles in this series

Addendum 20.10.2020: “Feature maps” or “response maps”?

Whilst reading some more books on AI I got the impression that some authors associate the term “feature map” with the full collection of all maps related to a Conv layer. A single map then is called a “response map” or an “output channel” of the feature map. However, as most books (e.g. the book of A. Geron on Machine Learning) call a singular map a “feature map” I cling to this usage throughout this series.

Samba 4, shares, wsdd and Windows 10 – how to list Linux Samba servers in the Win 10 Explorer

Posted on 24. May 2020 by eremo

These days I relatively often need to work with Windows 10 at home (home-office, corona virus, …). Normally, I isolate my own Win 10 instance in a VMware virtual machine on my Linux PC – and reduce any network connections of this VM to selected external servers. Under normal conditions all ports on the Linux host are closed for the virtual machine [VM]. But on a few temporary occasions I want to the Win 10 system to access a specific Samba exchange directory on a KVM virtualized Linux instance on the same host.

Off topic: You see that I never present directories of my Linux host directly to a Win 10 guest via Samba. Instead I transfer files via an exchange directory on an intermediate VM whose Samba service is configured to disallow access of the Win system on shares presented to the host. A primitive, but effective form of separation. The only inconvenient consequence is that synchronization becomes a two-fold process on the host and the Linux VM. But we have Linux tools for this, so the effort is limited. )

Of course we want to use the SMB protocol in a modern version, i.e. version 3.x (SMB3), over TCP/IP for this purpose (port 445). In addition we need some mechanism to detect and browse SMB servers on the Windows system. In the old days NetBIOS was used for the latter. On the Linux side we had the nmbd-daemon for it – and we could set up a special Samba server as a WINS server.

During the last year Microsoft has – via updates and new builds of Windows 10 – followed a consistent politics of deactivating the use of SMB V1.0 systematically. This, however, led to problems – not only between Windows PCs, but also between Win 10 instances and Samba 4 servers. This article addresses one of these problems: the missing list of available Samba servers in the Windows Explorer.

There are many contributions on the Internet describing this problem and some even say that you only can solve it by restoring SMB V1 capabilities in Win 10 again. In this article I want to recommend two different solutions:

Ignore the problem of Samba server detection and use your Samba shares on Win 10 with the SMB3 protocol as network drives.
If you absolutely want to see and list your Samba servers in the Windows Explorer of a Win 10 client, use the “Web-Service-Discovery” service via a WSDD-daemon provided by a Python script of Steffen Christgau.

I myself got on the right track of solving the named problem by an article of a guy called “Stilez”. His article is the first one listed under the section “Links” below. I recommend strongly to read it; it is Stilez who deserves all credit in pointing out both the problem and the solution. I just applied his insight to my own situation with virtualized Samba servers based on Opensuse Leap 15.1.

SMB V1.0 should be avoided – but NetBIOS needs it to exchange information about SMB servers

SMB, especially version SMB V1.0, is well known for security problems. Even MS has understood this – especially after the Wannacry disaster. See e.g. the links in the section “Links” => “Warnings of SMBV1” at the end of this article. MS has deactivated SMB V1 in the background via some updates of Win 8 and Win 10.

One of the resulting problem is that we do not see Samba servers in the Windows Explorer of a Win 10 system any longer. In the section “Network” of the Windows Explorer you normally should see a list of servers which are members of a Workgroup and offer shares.

Two years ago we would use NetBIOS’s discovery protocol and a WINS server to get this information. Unfortunately, the NetBIOS service detection ability depends on SMB1 features. The stupid thing is that we for a long while now had and have a relatively secure SMB2/3, but NetBIOS discovery only worked with SMB V1 enabled on the Windows client. Deactivating SMB V1 means deactivating NetBIOS at the
same time – and if you watch your Firewall logs for incoming packets from the Win 10 clients you will notice that exactly such a thing happened on Win 10 clients.

This actually means that you can have a full featured Samba/NetBIOS setup on the Linux side, that you may have opened the right ports on the firewalls for your Samba/WINS server and client systems, but that you will nevertheless not get any list of available Samba servers in Win 10’s Explorer. 🙁

Having understood this leads to the key question for our problem:

By what did MS replace the detection features of NetBIOS in combination with SMB-services?

Settings on the MS Win side – which alone will not help

When you google a bit you may find many hints regarding settings by which you activate network “discovery” functionalities via two Windows services. See

https://www.wintips.org/fix-windows-10-network-computers-not-showing/
https://winaero.com/blog/network-computers-not-visible-windows-10-version-1803/

You can follow these recommendations. If you want to see your own PC and other Windows systems in the Explorer’s list of network resources you must have activated them (see below). However, in my Win 10 client the recommended settings were already activated – with the exception of SMB V1, which I did and do not wish to reactivate again. The “discovery” settings may help you with other older Windows systems, but they do not enable a listing of Samba 4 servers without additional measures on Win 10.

There is another category of hints which in my opinion are contra-productive regarding security. See https://devanswers.co/network-error-problem-windows-cannot-access-hostname-samba/
Why activate an insecure setting? Especially, as such a setting does not really help us with our special problem? 🙁

A last set of hints concerns the settings on the Samba server, itself. I find it especially nice when such recommendations come from Microsoft :-). See: http://woshub.com/cannot-access-smb-network-shares-windows-10-1709/

[global]
server min protocol = SMB2_10
client max protocol = SMB3
client min protocol = SMB2_10
encrypt passwords = true
restrict anonymous = 2

Thanks to MS we now understand that we should not use SMB V1 …. But, actually, these hints are again insufficient regarding the Explorer problem …

What you could do – but should NOT do

Once you have understood that NetBIOS and SMB V1 still have an intimate relation (at least on a Windows systems) you may get the idea that there might exist some option to reactivate SMBV1 again on the Win 10 system. This is indeed possible. See here:
https://community.nethserver.org/t/windows-10-not-showing-servers-shares-in-network-browser/14263/4
https://www.wintips.org/fix-windows-10-network-computers-not-showing/

If you follow the advice of the authors and in addition re-open the standard ports for NetBIOS (UDP) 137, 138, (TCP) 139 on your firewalls between the Win 10 machine and your Samba servers you will – almost at once – get up the list of your accessible Samba servers in the Network section of the Win 10 Explorer. (Maybe you have to restart the smb and nmb services on your Linux machines).

But: You should not do this! SMB V1 should definitely become history!

Fortunately, a re-activation of SMB V1
on a Win 10 system is NOT required to mount Samba shares and it is neither required to get a list of available Samba servers in the Win 10 Explorer.

What you should do: Win 10 service settings

There are two service settings which are required to see other servers (and your own Win10 PC itself) in the list of network hosts presented by the Windows explorer:
Start services.msc ( press the Windows key + R => Enter “services.msc” in the dialog. Or: start services.msc it via the Control Panel => System and Security => Services)

Look for “Function Discovery Provider Host” => Set : Startup Type => Automatic
Look for “Function Discovery Resource Publication” => Set : Startup Type => Automatic (Delayed Start) !!

I noticed that on my VMware Win 10 guests the second setting appeared to be crucial to get the Win 10 PC itself listed among the network servers.

What you should do: Use the SMBV3 protocol!

As you as a Linux user meanwhile have probably replaced all your virtualized Win 7 guests, you should use the following settings in the [global] section of the configuration file “/etc/samba/smb.conf” of your Samba servers:

[global]
…
“protocol = SMB3”.
…

This is what Win 10 supports; you need SMB2_10 with some builds of Win 8 (???), only. Remember also that port 445 must be open on a firewall between the Win 10 client and your Samba server.

For Linux requirements to use SMB3 see
https://wiki.samba.org: SMB3 kernel status
For “SMB Direct” (RDMA) you normally need a kernel version > 4.16. On Opensuse Leap 15.1 most of the required kernel features have been backported. In Win 10 SMB Direct is normally activated; you find it in the “Window-Features” settings (https://www.windowscentral.com/how-manage-optional-features-windows-10)

Not seeing Samba servers in the Explorer does not mean that mounting a Samba share as a network drive does not work

Not seeing the Samba servers in the Win 10 Explorer – because the NetBIOS detection is defunct – does not mean that you cannot work with a Samba share on a Win 10 system. You can just “mount” it on Windows as a “network drive“:

Open a Windows Explorer, choose “This PC” on the left side, then click “Map network drive” in the upper area of the window and follow the instructions:
You choose a free drive letter and provide the Samba server name and its share in the usual MS form as “\\SERVERNAME\SHARE”.
Afterwards, you must activate the option “Connect using different credentials” in the dialog on the Win 10 side, if your Win 10 user for security reasons has a different UID and Password on the Samba server than on Win 10. Needless to say that this is a setting I strongly recommend – and of course we do not allow any direct anonymous or guest access to our Samba server without credentials delivered from a Windows machine (at least not without any central authentication systems).
So, you eventually must provide a valid Samba user name on your Samba server and the password – and there you happily go and use your resources on the Samba share from your Win 10 client.

I assumed of course that you have allowed access from the Win 10 host and the user by respective settings of “hosts allow” and “valid users” for the share in your Samba configuration.
Note: You need not mark the option for reconnecting the share in the Windows dialog for network drives if you only use the Samba exchange shares temporarily.

On an Opensuse system this works perfectly with the protocol settings for SMB3 on the server. So, you can use your shares even without seeing the samba
server in the Explorer: You just have to know what your shares are named and on which Samba servers they are located. No problem for a Linux admin.

In my opinion this approach is the most secure one among all “peer to peer”-approaches which have to work without a central network wide authentication service. It only requires to open port 445 for the time of a Samba session to a specific Samba server. Otherwise you do not provide any information for free to the Win 10 system and its “users”. (Well, an open question is what MS really does with the provided Samba credentials. But that is another story ….)

What you should do: Use the WSDD service on your Samba server

If you allow for some information sharing between your virtualized Win 10 and other KVM based virtual Samba machines in your LAN – and are not afraid of Microsoft or Antivirus companies on the Windows system to collect respective information – then there is a working option to get a stable list of the available Samba servers in the Windows Explorer – without the use of SMB V1.0.

Windows 10 implements web service detection via multiple mechanisms; among them: Multicast messages over ports 3702 (UDP), TCP 5357 and 1900 (UDP). For a detection of Samba services you “only” need ports 3702 (UDP) and 5357 (TCP). The general service detection port 1900 can remain closed in the firewalls between your Win 10 instances and your Linux world for our specific purpose. See
https://www.speedguide.net/port.php?port=5357
https://www.speedguide.net/port.php?port=3702
https://techcommunity.microsoft.com/t5/ask-the-performance-team/ws2008-the-wsd-port-monitor/ba-p/372760
https://en.wikipedia.org/wiki/Simple Service Discovery Protocol

The mechanism using ports 3702 and 5351 is called “Web Service Discovery” and was introduced by MS to cover the detection of printers and other devices in networks. In combination with SMB2 and SMB3 it is the preferred service to detect Samba services.

OK, do we have something like a counter-part available on a Linux system? Obviously, such a service is not (yet?) included in Samba 4 – at least not in the 4.9 version on my system with Opensuse Leap 15.1. The fact that WSD is not (yet?) a part of Samba may have some good reasons. See link.
One can understand the reservations and hesitation to include it, as WSD also serves other purposes than just the detection of SMB services.

Fortunately, a guy named Steffen Christgau, has written an (interesting) Python 3 script, which offers you the basic WSD functionality. See https://github.com/christgau/wsdd.

You can use the script in form of a daemon process on a Linux system – hence we speak of WSDD.

Using YaST I quickly found out that a WSDD RPM package is actually included in my “Opensuse Leap 15.1 Update” repository. People with other Linux distros may download the present WSDD version from GitHub.

On Opensuse it comes with an associated systemd service-file which you find in the directory “/usr/lib/systemd/system”.

[Unit]
Description=Web Services Dynamic Discovery host daemon
After=network-online.target
Wants=network-online.target

[Service]
Type=simple
AmbientCapabilities=CAP_SYS_CHROOT
PermissionsStartOnly=true
Environment= WSDD_ARGS=-p
ExecStartPre=/usr/lib/wsdd/wsdd-init.sh
EnvironmentFile=-/run/sysconfig/wsdd
ExecStart=/usr/sbin/wsdd --shortlog -c /run/wsdd $WSDD_ARGS
ExecStartPost=/usr/bin/rm /run/
sysconfig/wsdd
User=wsdd
Group=wsdd

[Install]
WantedBy=multi-user.target

Reading the documentation you find out that the daemon runs chrooted – which is a reasonable security measure.
Opensuse even provides an elementary configuration file in “/etc/sysconfig/wsdd“.

I used the parameter

WSDD_WORKGROUP=”MYWORKGROUP”

there to announce the right Workgroup for my (virtualized) Samba server.

So, I had everything ready to start WSDD by “rcwsdd start” (or by “systemctl start wsdd.service”) on my Samba server.

On the local firewall of the SMB server I opened

port 445 (TCP) for SMB(3) In/Out for the server and from/to the Win-10-Client,
port 3702 (UDP) for incoming packets to the server and outgoing packets from the server to the Multicast address 239.255.255.250,
port 5357 (TCP) In/Out for the server and from/to the Win 10 client.

And: I closed all NetBIOS ports (UDP 137, 138 / TCP 139) and eventually stopped the “nmbd”-service on the Samba server! (UDP 137, 138 / TCP 139)

Within a second or so, my Samba 4 server appeared in the Windows 10 Explorer!

Further hints:
As the 3702 port is used with the UDP protocol it should be regarded as potentially dangerous. See: https://blogs.akamai.com/sitr/2019/09/new-ddos-vector-observed-in-the-wild-wsd-attacks-hitting-35gbps.html
The port 1900 which appeared in the firewall logs does not seem to be important. I blocked it.

So far, so good. However, when I refreshed the list in the Win 10 Explorer my SAMBA server disappeared again. 🙁

What you should do: Take special care about the network interface to which the WSDD service gets attached to

It took me a while to find out that the origin of the last problem had to do with the fact that my virtualized server and my Win 10 client both had multiple network interfaces on virtualized bridges. There are no loops in the configuration, but it occurred that multiple broadcasts packets arrive via different paths at the Samba server and were answered – and thus multiple return messages appeared at the Win 10 client during a refresh – which Win 10 did not like (see the discussion in the following link.
https://github.com/christgau/wsdd/issues/8

As soon as I restricted the answer of the Samba server to exactly one of the interfaces on my virtual bridge via the the parameter “WSDD_INTERFACES” in the “/etc/sysconfig/wsdd”-configuration file everything went fine. Refreshes now lead to an immediate update including the Samba server.

So, be a little careful, when you have some complicated bridge structures associated with your virtualized VMware or KVM guests. The WSDD service should be limited to exactly one interface of the Samba server.

Note: As we do not need NetBIOS any longer – block ports 137, 138 (UDP) and 139 (TCP) in your firewalls! It will make you feel better instantaneously.

Conclusion

The “end” of SMB V1 on Win 10 is a reasonable step. However, it undermines the visibility of Samba servers in the Windows Explorers. The reason is that NetBIOS requires SMB1.0 features on Windows. NetBIOS is/was therefore consistently deactivated on Win 10, too. The service detection on the network is replaced by the WSD service which was originally introduced for printer detection (and possibly other devices). Activating it on the Win 10 system may help with the detection of other Windows (8 and 10) systems on the network, but not with Samba 4 servers. Samba servers presently only serve NetBIOS requests of Win clients
to allow for server and share detection. Therefore, without additional measures, they are not displayed in the Windows Explorer of a regular Win 10 client.

This does, however, not restrict the usage of Samba shares on the Win 10 client via the SMB3 protocol. They can be used as “network drives” – just as before. Not distributing name and device information on a network has its advantages regarding security.

If you absolutely must see your Samba servers in the Win 10 Explorer install and configure the WSDD package of Steffen Christgau. You can use it as a systemd service. You should restrict the interfaces WSDD gets attached to – especially if your Samba servers are attached to virtual network bridges (Linux bridges or VMware bridges).

So:

Disable SMBV1 in Windows 10 if an update has not yet done it for you!
Set the protocol in the Samba servers to SMBV3!
Try to work with “networks drives” on your Win 10 guests, only!
Install, configure and use WSDD, if you really need to see your Samba servers in the Windows Explorer.
Open the port 445 (TCP, IN/OUT between the Win 10 client and the server), 3072 (UDP, OUT from the server and the Win 10 client to 239.255.255.250, IN to the server from the Win 10 client / IN to the Win 10 client from the server; rules details depending on the firewall location), port 5357 (TCP; In/OUT between the Samba server and the Win 10 client) on your firewalls between the Samba server and the Win 10 system.
Close the NetBIOS ports in your firewalls!
You should also take care of stopping multicast messages leaving perimeter firewalls; normally packets to multicast addresses should not be routed, but blocking them explicitly for certain interfaces is no harm, either.

Of course you must repeat the WSDD and firewall setup for all your Samba servers. But as a Linux admin you have your tools for distributing common configuration files or copying virtualization setups.

Links

The real story
!!!! https://www.ixsystems.com/community/resources/how-to-kill-off-smb1-netbios-wins-and-still-have-windows-network-neighbourhood-better-than-ever.106/ !!!

https://forums.linuxmint.com/viewtopic.php?p=1799875

https://devanswers.co/discover-ubuntu-machines-samba-shares-windows-10-network/

https://bugs.launchpad.net/ubuntu/ source/ samba/ +bug/ 1831441

https://forums.opensuse.org/ showthread.php/ 540083-Samba-Network-Device-Type-for-Windows-10

https://kofler.info/zugriff-auf-netzwerkverzeichnisse-mit-nautilus/

WSDD and its problems
https://github.com/christgau/wsdd
https://github.com/christgau/wsdd/issues/8
https://forums.opensuse.org/ showthread.php/ 540083-Samba-Network-Device-Type-for-Windows-10

Warnings of SMB V1
https://docs.microsoft.com/de-de/windows-server/storage/file-server/troubleshoot/detect-enable-and-disable-smbv1-v2-v3
https://blog.malwarebytes.com/101/2018/12/how-threat-actors-are-using-smb-vulnerabilities/
https://securityboulevard.com/2018/12/whats-the-problem-with-smb-1-and-should-you-worry-about-smb-2-and-3/
https://techcommunity.microsoft.com/t5/storage-at-microsoft/stop-using-smb1/ba-p/425858
https://www.cubespotter.de/cubespotter/wannacry-nsa-exploits-und-das-maerchen-von-smbv1/

Problems with Win 10 and shares
https://social.technet.microsoft.com/ Forums/ en-US: cannot-connect-to-cifs-smb-samba-network-shares-amp-shared-folders-in-windows-10-after-update?forum=win10itpronetworking

RDMA and SMB Direct
https://searchstorage.techtarget.com/ definition/ Remote-Direct-Memory-Access

Other settings in the SMB/Samba environment of minor relevance
http://woshub.com/cannot-access-smb-network-shares-windows-10-1709/
https://superuser.com/questions/1466968/unable-to-connect-to-a-linux-samba-server-via-hostname-on-windows-10
https://superuser.com/questions/1522896/windows-10-cannot-connect-to-linux-samba-shares-except-from-smb1-cifs
https://www.reddit.com/ r/ techsupport/ comments/ 3yevip/ windows 10 cant see samba shares/
https://devanswers.co/network-error-problem-windows-cannot-access-hostname-samba/

MLP, Numpy, TF2 – performance issues – Step III – a correction to BW propagation

Posted on 15. May 2020 by eremo

In the last articles of this series

MLP, Numpy, TF2 – performance issues – Step II – bias neurons, F- or C- contiguous arrays and performance
MLP, Numpy, TF2 – performance issues – Step I – float32, reduction of back propagation

we looked at the FW-propagation of the MLP code which I discussed in another article series. We found that the treatment of bias neurons in the input layer was technically inefficient due to a collision of C- and F-contiguous arrays. By circumventing the problem we could accelerate the FW-propagation of big batches (as the complete training or test data set) by more than a factor of 2.

In this article I want to turn to the BW propagation and do some analysis regarding CPU consumption there. We will find a simple (and stupid) calculation step there which we shall replace. This will give us another 15% to 22% performance improvement in comparison to what we have reached in the last article for MNIST data:

9.6 secs for 35 epochs and a batch-size of 500
and 8.7 secs for a batch-size of 20000.

Present CPU time relation between the FW- and the BW-propagation

The central training of mini-batches is performed by the method “_handle_mini_batch()”.

#
    ''' -- Method to deal with a batch -- '''
    def _handle_mini_batch (self, num_batch = 0, num_epoch = 0, b_print_y_vals = False, b_print = False, b_keep_bw_matrices = True):
        ''' .... '''
        # Layer-related lists to be filled with 2-dim Numpy matrices during FW propagation
        # ********************************************************************************
        li_Z_in_layer  = [None] * self._n_total_layers # List of matrices with z-input values for each layer; filled during FW-propagation
        li_A_out_layer = li_Z_in_layer.copy()          # List of matrices with results of activation/output-functions for each layer; filled during FW-propagation
        li_delta_out   = li_Z_in_layer.copy()          # Matrix with out_delta-values at the outermost layer 
        li_delta_layer = li_Z_in_layer.copy()          # List of the matrices for the BW propagated delta values 
        li_D_layer     = li_Z_in_layer.copy()          # List of the derivative matrices D containing partial derivatives of the activation/ouput functions 
        li_grad_layer  = li_Z_in_layer.copy()          # List of the matrices with gradient values for weight corrections
        
        # Major steps for the mini-batch during one epoch iteration 
        # **********************************************************
        
        #ts=time.perf_counter()
        # Step 0: List of indices for data records in the present mini-batch
        # ******
        ay_idx_batch = self._ay_mini_batches[num_batch]
        
        # Step 1: Special preparation of the Z-input to the MLP's input Layer L0
        # ******
        # ts=time.perf_counter()
        # slicing 
        li_Z_in_layer[0] = self._X_train[ay_idx_batch] # numpy arrays can be indexed by an array of integers
        
        # transposition 
        #~~~~~~~~~~~~~~
        li_Z_in_layer[0] = li_Z_in_layer[0].T
        #te=time.perf_counter(); t_batch = te - ts;
        #print("\nti - transposed inputbatch =", t_batch)
        
        # Step 2: Call forward propagation method for the present mini-batch of training records
        # *******
n        #tsa = time.perf_counter() 
        self._fw_propagation(li_Z_in = li_Z_in_layer, li_A_out = li_A_out_layer) 
        #tea = time.perf_counter(); ta = tea - tsa;  print("ta - FW-propagation", "%10.8f"%ta)
        
        # Step 3: Cost calculation for the mini-batch 
        # ********
        #tsb = time.perf_counter() 
        ay_y_enc = self._ay_onehot[:, ay_idx_batch]
        ay_ANN_out = li_A_out_layer[self._n_total_layers-1]
        total_costs_batch, rel_reg_contrib = self._calculate_loss_for_batch(ay_y_enc, ay_ANN_out, b_print = False)
        # we add the present cost value to the numpy array 
        self._ay_costs[num_epoch, num_batch]            = total_costs_batch
        self._ay_reg_cost_contrib[num_epoch, num_batch] = rel_reg_contrib
        #teb = time.perf_counter(); tb = teb - tsb; print("tb - cost calculation", "%10.8f"%tb)
        
        
        # Step 4: Avg-error for later plotting 
        # ********
        #tsc = time.perf_counter() 
        # mean "error" values - averaged over all nodes at outermost layer and all data sets of a mini-batch 
        ay_theta_out = ay_y_enc - ay_ANN_out
        ay_theta_avg = np.average(np.abs(ay_theta_out)) 
        self._ay_theta[num_epoch, num_batch] = ay_theta_avg 
        #tec = time.perf_counter(); tc = tec - tsc; print("tc - error", "%10.8f"%tc)
        
        
        # Step 5: Perform gradient calculation via back propagation of errors
        # ******* 
        #tsd = time.perf_counter() 
        self._bw_propagation( ay_y_enc = ay_y_enc, 
                              li_Z_in = li_Z_in_layer, 
                              li_A_out = li_A_out_layer, 
                              li_delta_out = li_delta_out, 
                              li_delta = li_delta_layer,
                              li_D = li_D_layer, 
                              li_grad = li_grad_layer, 
                              b_print = b_print,
                              b_internal_timing = False 
                              ) 
        #ted = time.perf_counter(); td = ted - tsd; print("td - BW propagation", "%10.8f"%td)
        
        # Step 7: Adjustment of weights  
        # *******        
        #tsf = time.perf_counter() 
        rg_layer=range(0, self._n_total_layers -1)
        for N in rg_layer:
            delta_w_N = self._learn_rate * li_grad_layer[N]
            self._li_w[N] -= ( delta_w_N + (self._mom_rate * self._li_mom[N]) )
            
            # save momentum
            self._li_mom[N] = delta_w_N
        #tef = time.perf_counter(); tf = tef - tsf; print("tf - weight correction", "%10.8f"%tf)
        
        return None

I took some time measurements there:

ti - transposed inputbatch = 0.0001785
ta - FW-propagation 0.00080975
tb - cost calculation 0.00030705
tc - error 0.00016182
td - BW propagation 0.00112558
tf - weight correction 0.00020079

ti - transposed inputbatch = 0.00018144
ta - FW-propagation 0.00082022
tb - cost calculation 0.00031284
tc - error 0.00016652
td - BW propagation 0.00106464
tf - weight correction 0.00019576

You see that the FW-propagation is a bit faster than the BW-propagation. This is a bit strange as the FW-propagation is dominated meanwhile by a really expensive operation which we cannot accelerate (without choosing a new activation function): The calculation of the sigmoid value for the inputs at layer L1.

So let us look into the BW-propagation; the code for it is momentarily:

    ''' -- Method to handle error BW propagation for a mini-batch --'''
    def _bw_propagation(self, 
                        ay_y_enc, li_Z_in, li_A_out, 
                        li_delta_out, li_delta, li_D, li_
grad, 
                        b_print = True, b_internal_timing = False):
        
        # List initialization: All parameter lists or arrays are filled or to be filled by layer operations 
        # Note: the lists li_Z_in, li_A_out were already filled by _fw_propagation() for the present batch 
        
        # Initiate BW propagation - provide delta-matrices for outermost layer
        # *********************** 
        tsa = time.perf_counter() 
        # Input Z at outermost layer E  (4 layers -> layer 3)
        ay_Z_E = li_Z_in[self._n_total_layers-1]
        # Output A at outermost layer E (was calculated by output function)
        ay_A_E = li_A_out[self._n_total_layers-1]
        
        # Calculate D-matrix (derivative of output function) at outmost the layer - presently only D_sigmoid 
        ay_D_E = self._calculate_D_E(ay_Z_E=ay_Z_E, b_print=b_print )
        #ay_D_E = ay_A_E * (1.0 - ay_A_E)

        # Get the 2 delta matrices for the outermost layer (only layer E has 2 delta-matrices)
        ay_delta_E, ay_delta_out_E = self._calculate_delta_E(ay_y_enc=ay_y_enc, ay_A_E=ay_A_E, ay_D_E=ay_D_E, b_print=b_print) 
        
        # add the matrices to their lists ; li_delta_out gets only one element 
        idxE = self._n_total_layers - 1
        li_delta_out[idxE] = ay_delta_out_E # this happens only once
        li_delta[idxE]     = ay_delta_E
        li_D[idxE]         = ay_D_E
        li_grad[idxE]      = None    # On the outermost layer there is no gradient ! 
        
        tea = time.perf_counter(); ta = tea - tsa; print("\nta-bp", "%10.8f"%ta)
        
        # Loop over all layers in reverse direction 
        # ******************************************
        # index range of target layers N in BW direction (starting with E-1 => 4 layers -> layer 2))
        range_N_bw_layer = reversed(range(0, self._n_total_layers-1))   # must be -1 as the last element is not taken 
        
        # loop over layers 
        tsb = time.perf_counter() 
        for N in range_N_bw_layer:
            
            # Back Propagation operations between layers N+1 and N 
            # *******************************************************
            # this method handles the special treatment of bias nodes in Z_in, too
            tsib = time.perf_counter() 
            ay_delta_N, ay_D_N, ay_grad_N = self._bw_prop_Np1_to_N( N=N, li_Z_in=li_Z_in, li_A_out=li_A_out, li_delta=li_delta, b_print=False )
            teib = time.perf_counter(); tib = teib - tsib; print("N = ", N, " tib-bp", "%10.8f"%tib)
            
            # add matrices to their lists 
            #tsic = time.perf_counter() 
            li_delta[N] = ay_delta_N
            li_D[N]     = ay_D_N
            li_grad[N]= ay_grad_N
            #teic = time.perf_counter(); tic = teic - tsic; print("\nN = ", N, " tic = ", "%10.8f"%tic)
        teb = time.perf_counter(); tb = teb - tsb; print("tb-bp", "%10.8f"%tb)
       
        return

Typical timing results are:

ta-bp 0.00007112
N =  2  tib-bp 0.00025399
N =  1  tib-bp 0.00051683
N =  0  tib-bp 0.00035941
tb-bp 0.00126436

ta-bp 0.00007492
N =  2  tib-bp 0.00027644
N =  1  tib-bp 0.00090043
N =  0  tib-bp 0.00036728
tb-bp 0.00168378

We see that the CPU consumption of “_bw_prop_Np1_to_N()” should be analyzed in detail. It is relatively time consuming at every layer, but especially at layer L1. (The list adds are insignificant.)
What does this method presently look like?

    ''' -- Method to calculate the BW-propagated delta-matrix and the gradient matrix to/for layer N '''
    def _bw_prop_Np1_to_N(self, N, li_Z_in, li_A_out, li_delta, b_print=False):
        '''
        BW-
error-propagation between layer N+1 and N 
        Version 1.5 - partially accelerated 

        Inputs: 
            li_Z_in:  List of input Z-matrices on all layers - values were calculated during FW-propagation
            li_A_out: List of output A-matrices - values were calculated during FW-propagation
            li_delta: List of delta-matrices - values for outermost ölayer E to layer N+1 should exist 
        
        Returns: 
            ay_delta_N - delta-matrix of layer N (required in subsequent steps)
            ay_D_N     - derivative matrix for the activation function on layer N 
            ay_grad_N  - matrix with gradient elements of the cost fnction with respect to the weights on layer N 
        '''
        
        # Prepare required quantities - and add bias neuron to ay_Z_in 
        # ****************************
        
        # Weight matrix meddling between layers N and N+1 
        ay_W_N = self._li_w[N]

        # delta-matrix of layer N+1
        ay_delta_Np1 = li_delta[N+1]

        # fetch output value saved during FW propagation 
        ay_A_N = li_A_out[N]

        # Optimization V1.5 ! 
        if N > 0: 
            
            #ts=time.perf_counter()
            ay_Z_N = li_Z_in[N]
            # !!! Add intermediate row (for bias) to Z_N !!!
            ay_Z_N = self._add_bias_neuron_to_layer(ay_Z_N, 'row')
            #te=time.perf_counter(); t1 = te - ts; print("\nBW t1 = ", t1, " N = ", N) 
        
            # Derivative matrix for the activation function (with extra bias node row)
            # ********************
            #    can only be calculated now as we need the z-values
            #ts=time.perf_counter()
            ay_D_N = self._calculate_D_N(ay_Z_N)
            #te=time.perf_counter(); t2 = te - ts; print("\nBW t2 = ", t2, " N = ", N) 
            
            # Propagate delta
            # **************

            # intermediate delta 
            # ~~~~~~~~~~~~~~~~~~
            #ts=time.perf_counter()
            ay_delta_w_N = ay_W_N.T.dot(ay_delta_Np1)
            #te=time.perf_counter(); t3 = te - ts; print("\nBW t3 = ", t3) 
            
            # final delta 
            # ~~~~~~~~~~~
            #ts=time.perf_counter()
            ay_delta_N = ay_delta_w_N * ay_D_N
            
            # Orig reduce dimension again
            # **************************** 
            ay_delta_N = ay_delta_N[1:, :]
            #te=time.perf_counter(); t4 = te - ts; print("\nBW t4 = ", t4) 
            
        else: 
            ay_delta_N = None
            ay_D_N = None
        
        # Calculate gradient
        # ********************
        #ts=time.perf_counter()
        ay_grad_N = np.dot(ay_delta_Np1, ay_A_N.T)
        #te=time.perf_counter(); t5 = te - ts; print("\nBW t5 = ", t5) 
        
        # regularize gradient (!!!! without adding bias nodes in the L1, L2 sums) 
        #ts=time.perf_counter()
        if self._lambda2_reg > 0: 
            ay_grad_N[:, 1:] += self._li_w[N][:, 1:] * self._lambda2_reg 
        if self._lambda1_reg > 0: 
            ay_grad_N[:, 1:] += np.sign(self._li_w[N][:, 1:]) * self._lambda1_reg 
        #te=time.perf_counter(); t6 = te - ts; print("\nBW t6 = ", t6) 
        
        return ay_delta_N, ay_D_N, ay_grad_N

Timing data for a batch-size of 500 are:

N =  2
BW t1 =  0.0001169009999557602  N =  2
BW t2 =  0.00035331499998392246  N =  2
BW t3 =  0.00018078099992635543
BW t4 =  0.00010234199999104021
BW t5 =  9.928200006470433e-05
BW t6 =  2.4267000071631628e-05
N =  2  tib-bp 0.00124414

N =  1
BW t1 =  0.0004323499999827618  N =  1
BW t2 =  0.
000781415999881574  N =  1
BW t3 =  4.2077999978573644e-05
BW t4 =  0.00022921000004316738
BW t5 =  9.376399998473062e-05
BW t6 =  0.00012183700005152787
N =  1  tib-bp 0.00216281

N =  0
BW t5 =  0.0004289769999559212
BW t6 =  0.00015404999999191205
N =  0  tib-bp 0.00075249
....
N =  2
BW t1 =  0.00012802800006284087  N =  2
BW t2 =  0.00034988200013685855  N =  2
BW t3 =  0.0001854429999639251
BW t4 =  0.00010359299994888715
BW t5 =  0.00010210400000687514
BW t6 =  2.4010999823076418e-05
N =  2  tib-bp 0.00125854

N =  1
BW t1 =  0.0004407169999467442  N =  1
BW t2 =  0.0007845899999665562  N =  1
BW t3 =  0.00025684100000944454
BW t4 =  0.00012409999999363208
BW t5 =  0.00010345399982725212
BW t6 =  0.00012994100006835652
N =  1  tib-bp 0.00221321

N =  0
BW t5 =  0.00044504700008474174
BW t6 =  0.00016473000005134963
N =  0  tib-bp 0.00071442

....
N =  2
BW t1 =  0.000292730999944979  N =  2
BW t2 =  0.001102525000078458  N =  2
BW t3 =  2.9429999813146424e-05
BW t4 =  8.547999868824263e-06
BW t5 =  3.554099998837046e-05
BW t6 =  2.5041999833774753e-05
N =  2  tib-bp 0.00178565

N =  1
BW t1 =  3.143399999316898e-05  N =  1
BW t2 =  0.0006720640001276479  N =  1
BW t3 =  5.4785999964224175e-05
BW t4 =  9.756200006449944e-05
BW t5 =  0.0001605449999715347
BW t6 =  1.8391000139672542e-05
N =  1  tib-bp 0.00147566

N =  0
BW t5 =  0.0003641810001226986
BW t6 =  6.338999992294703e-05
N =  0  tib-bp 0.00046542

It seems that we should care about t1, t2, t3 for hidden layers and maybe about t5 at layers L1/L0.

However, for a batch-size of 15000 things look a bit different:

N =  2
BW t1 =  0.0005776280000304723  N =  2
BW t2 =  0.004995969999981753  N =  2
BW t3 =  0.0003165199999557444
BW t4 =  0.0005244750000201748
BW t5 =  0.000518499999998312
BW t6 =  2.2458999978880456e-05
N =  2  tib-bp 0.00736144

N =  1
BW t1 =  0.0010120430000029046  N =  1
BW t2 =  0.010797029000002567  N =  1
BW t3 =  0.0005006920000028003
BW t4 =  0.0008704929999794331
BW t5 =  0.0010805200000163495
BW t6 =  3.0326000000968634e-05
N =  1  tib-bp 0.01463436

N =  0
BW t5 =  0.006987539000022025
BW t6 =  0.00023552499999368592
N =  0  tib-bp 0.00730959


N =  2
BW t1 =  0.0006299790000525718  N =  2
BW t2 =  0.005081416999985322  N =  2
BW t3 =  0.00018547400003399162
BW t4 =  0.0005970070000103078
BW t5 =  0.000564008000026206
BW t6 =  2.3311000006742688e-05
N =  2  tib-bp 0.00737899

N =  1
BW t1 =  0.0009376909999900818  N =  1
BW t2 =  0.010650266999959968  N =  1
BW t3 =  0.0005232729999988806
BW t4 =  0.0009100700000317374
BW t5 =  0.0011237720000281115
BW t6 =  0.00016643800000792908
N =  1  tib-bp 0.01466144

N =  0
BW t5 =  0.006987463000029948
BW t6 =  0.00023978600000873485
N =  0  tib-bp 0.00734308

For big batch-sizes “t2” dominates everything. It seems that we have found another code area which causes the trouble with big batch-sizes which we already observed before!

What operations do the different CPU times stand for?

To keep an overview without looking into the code again, I briefly summarize which operations cause which of the measured time differences:

“t1” – which contributes for small batch-sizes stands for adding a bias neuron to the input data Z_in at each layer.
“t2” – which is by far dominant for big batch sizes stands for calculating the derivative of the output/activation function (in our case of the sigmoid function) at the various layers.
“t3” – which contributes at
some layers stands for a dot()-matrix multiplication with the transposed weight-matrix,
“t4” – covers an element-wise matrix-multiplication,
“t5” – contributes at the BW-transition from layer L1 to L0 and covers the matrix multiplication there (including the full output matrix with the bias neurons at L0)

Use the output values calculated at each layer during FW-propagation!

Why does the calculation of the derivative of the sigmoid function take so much time? Answer: Because I coded it stupidly! Just look at it:

    ''' -- Method to calculate the matrix with the derivative values of the output function at outermost layer '''
    def _calculate_D_N(self, ay_Z_N, b_print= False):
        '''
        This method calculates and returns the D-matrix for the outermost layer
        The D matrix contains derivatives of the output function with respect to local input "z_j" at outermost nodes. 
        
        Returns
        ------
        ay_D_E:    Matrix with derivative values of the output function 
                   with respect to local z_j valus at the nodes of the outermost layer E
        Note: This is a 2-dim matrix over layer nodes and training samples of the mini-batch
        '''
        if self._my_out_func == 'sigmoid':
            ay_D_E = self._D_sigmoid(ay_Z = ay_Z_N)
        
        else:
            print("The derivative for output function " + self._my_out_func + " is not known yet!" )
            sys.exit()
        
        return ay_D_E

    ''' -- method for the derivative of the sigmoid function-- '''
    def _D_sigmoid(self, ay_Z):
        ''' 
        Derivative of sigmoid function with respect to Z-values 
        - works via expit element-wise on matrices
        Input:  Z - Matrix with Input values for the activation function Phi() = sigmoid() 
        Output: D - Matrix with derivative values 
        '''
        S_Z = self._sigmoid(ay_Z)
        return S_Z * (1.0 - S_Z)

We first call an intermediate function which then directs us to the right function for a chosen activation function. Well meant: So far, we use only the sigmoid function, but it could e.g. also be the relu() or tanh()-function. So, we did what we did for the sake of generalization. But we did it badly because of two reasons:

We did not keep up a function call pattern which we introduced in the FW-propagation.
The calculation of the derivative is inefficient.

The first point is a minor one: During FW-propagation we called the right (!) activation function, i.e. the one we choose by input parameters to our ANN-object, by an indirect call. Why not do it the same way here? We would avoid an intermediate function call and keep up a pattern. Actually, we prepared the necessary definitions already in the __init__()-function.

The second point is relevant for performance: The derivative function produces the correct results for a given “ay_Z”, but this is totally inefficient in our BW-situation. The code repeats a really expensive operation which we have already performed during FW-propagation: calling sigmoid(ay_Z) to get “A_out”-values per layer then. We even put the A_out-values [=sigmoid(ay_Z_in)] per layer and batch (!) with some foresight into a list in “li_A_out[]” at that point of the code (see the FW-propagation code discussed in the last article).

So, of course, we should use these “A_out”-values now in the BW-steps! No further comment …. you see what we need to do.

Hint: Actually, also other activation functions “act(Z)” like e.g. the “tanh()”-function have derivatives which depend on on “A=act(
Z)”, only. So, we should provide Z and A via an interface to the derivative function and let the respective functions take what it needs.
But, my insight into my own dumbness gets worse.

Eliminate the bias neuron operation!

Why did we need a bias-neuron operation? Answer: We do not need it! It was only introduced due to insufficient cleverness. In the article

A simple Python program for an ANN to cover the MNIST dataset – VII – EBP related topics and obstacles

I have already indicated that we use the function for adding a row of bias-neurons again only to compensate one deficit: The matrix of the derivative values did not fit the shape of the weight matrix for the required element-wise operations. However, I also said: There probably is an alternative.

Well, let me make a long story short: The steps behind t1 up to t4 to calculate “ay_delta_N” for the present layer L_N (with N>=1) can be compressed into two relatively simple lines:

ay_delta_w_N = ay_W_N.T.dot(ay_delta_Np1)
ay_delta_N = ay_delta_w_N[1:,:] * ay_A_N[1:,:] * (1.0 – ay_A_N[1:,:]); ay_D_N = None;

No bias back and forth corrections! Instead we use simple slicing to compensate for our weight matrices with a shape covering an extra row of bias node output. No Z-based derivative calculation; no sigmoid(Z)-call. The last statement is only required to support the present output interface. Think it through in detail; the shortcut does not cause any harm.

Code change for tests

Before we bring the code into a new consolidated form with re-coded methods let us see what we gain by just changing the code to the two lines given above in terms of CPU time and performance. Our function “_bw_prop_Np1_to_N()” then gets reduced to the following lines:

    ''' -- Method to calculate the BW-propagated delta-matrix and the gradient matrix to/for layer N '''
    def _bw_prop_Np1_to_N(self, N, li_Z_in, li_A_out, li_delta, b_print=False):
        
        # Weight matrix meddling between layers N and N+1 
        ay_W_N = self._li_w[N]
        ay_delta_Np1 = li_delta[N+1]

        # fetch output value saved during FW propagation 
        ay_A_N = li_A_out[N]

        # Optimization from previous version  
        if N > 0: 
            #ts=time.perf_counter()
            ay_Z_N = li_Z_in[N]
            
            # Propagate delta
            # ~~~~~~~~~~~~~~~~~
            ay_delta_w_N = ay_W_N.T.dot(ay_delta_Np1)
            ay_delta_N = ay_delta_w_N[1:,:] * ay_A_N[1:,:] * (1.0 - ay_A_N[1:,:])
            ay_D_N = None; 
            
        else: 
            ay_delta_N = None
            ay_D_N = None
        
        # Calculate gradient
        # ********************
        ay_grad_N = np.dot(ay_delta_Np1, ay_A_N.T)
        
        if self._lambda2_reg > 0: 
            ay_grad_N[:, 1:] += self._li_w[N][:, 1:] * self._lambda2_reg 
        if self._lambda1_reg > 0: 
            ay_grad_N[:, 1:] += np.sign(self._li_w[N][:, 1:]) * self._lambda1_reg 
        
        return ay_delta_N, ay_D_N, ay_grad_N

Performance gain

What run times do we get with this setting? We perform our typical test runs over 35 epochs – but this time for two different batch-sizes:

Batch-size = 500

 
------------------
Starting epoch 35

Time_CPU for epoch 35 0.2169024469985743
Total CPU-time:  7.52385053600301

learning rate =  0.0009994051838157095

total costs of training set   =  -1.0
rel. reg. contrib. to total costs =  -1.0

total costs 
of last mini_batch   =  65.43618
rel. reg. contrib. to batch costs =  0.12302863

mean abs weight at L0 :  -10.0
mean abs weight at L1 :  -10.0
mean abs weight at L2 :  -10.0

avg total error of last mini_batch =  0.00758
presently batch averaged accuracy   =  0.99272

-------------------
Total training Time_CPU:  7.5257336139984545

Not bad! We became faster by around 2 secs compared to the results of the last article! This is close to an improvement of 20%.

But what about big batch sizes? Here is the result for a relatively big batch size:

Batch-size = 20000

------------------
Starting epoch 35

Time_CPU for epoch 35 0.2019189490019926
Total CPU-time:  6.716679593999288

learning rate =  9.994051838157101e-05

total costs of training set   =  -1.0
rel. reg. contrib. to total costs =  -1.0

total costs of last mini_batch   =  13028.141
rel. reg. contrib. to batch costs =  0.00021923862

mean abs weight at L0 :  -10.0
mean abs weight at L1 :  -10.0
mean abs weight at L2 :  -10.0

avg total error of last mini_batch =  0.04389
presently batch averaged accuracy   =  0.95602

-------------------
Total training Time_CPU:  6.716954112998792

Again an acceleration by roughly 2 secs – corresponding to an improvement of 22%!

In both cases I took the best result out of three runs.

Conclusion

Enough for today! We have done a major step with regard to performance optimization also in the BW-propagation. It remains to re-code the derivative calculation in form which uses indirect function calls to remain flexible. I shall give you the code in the next article.

We learned today is that we, of course, should reuse the results of the FW-propagation and that it is indeed a good investment to save the output data per layer in some Python list or other suitable structures during FW-propagation. We also saw again that a sufficiently efficient bias neuron treatment can be achieved by a more efficient solution than provisioned so far.

All in all we have meanwhile gained more than a factor of 6.5 in performance since we started with optimization. Our new standard values are 7.3 secs and 6.8 secs for 35 epochs on MNIST data and batch sizes of 500 and 20000, respectively.

We have reached the order of what Keras and TF2 can deliver on a CPU for big batch sizes. For small batch sizes we are already faster. This indicates that we have done no bad job so far …

In the next article we shall look a bit at the matrix operations and evaluate further optimization options.

Linux-Blog – Dr. Mönchmeyer / anracon

Notes about Linux, ML and some simple math …

Author Archives: eremo

A simple CNN for the MNIST datasets – I – CNN basics

Why CNNs?

Filters, kernels and feature maps

Filters and their related weights are the end products of the training and optimization of a CNN!

Nothing is really “moved” …

Activation functions

Tensors

Chained convolutions cover growing areas of the original image

Sub-sampling is necessary!

Resulting options for architectures

Further aspects

Questions and objectives for this article series

Conclusion

Literature

Further articles in this series

Addendum 20.10.2020: “Feature maps” or “response maps”?

Samba 4, shares, wsdd and Windows 10 – how to list Linux Samba servers in the Win 10 Explorer

SMB V1.0 should be avoided – but NetBIOS needs it to exchange information about SMB servers

Settings on the MS Win side – which alone will not help

What you could do – but should NOT do

What you should do: Win 10 service settings

What you should do: Use the SMBV3 protocol!

Not seeing Samba servers in the Explorer does not mean that mounting a Samba share as a network drive does not work

What you should do: Use the WSDD service on your Samba server

What you should do: Take special care about the network interface to which the WSDD service gets attached to

Conclusion

Links

MLP, Numpy, TF2 – performance issues – Step III – a correction to BW propagation

Present CPU time relation between the FW- and the BW-propagation

What operations do the different CPU times stand for?

Use the output values calculated at each layer during FW-propagation!

Eliminate the bias neuron operation!

Code change for tests

Performance gain

Conclusion