# Towards Better Input Masking for Convolutional Neural Networks

Sriram Balasubramanian
University of Maryland, College Park
<EMAIL_ADDRESS>Soheil Feizi
University of Maryland, College Park
<EMAIL_ADDRESS>
###### Abstract

The ability to remove features from the input of machine learning models is
very important to understand and interpret model predictions. However, this is
non-trivial for vision models since masking out parts of the input image and
replacing them with a baseline color like black or grey typically causes large
distribution shifts. Masking may even make the model focus on the masking
patterns for its prediction rather than the unmasked portions of the image. In
recent work, it has been shown that vision transformers are less affected by
such issues as one can simply drop the tokens corresponding to the masked
image portions. They are thus more easily interpretable using techniques like
LIME which rely on input perturbation. Using the same intuition, we devise a
masking technique for CNNs called layer masking, which simulates running the
CNN on only the unmasked input. We find that our method is (i) much less
disruptive to the model’s output and its intermediate activations, and (ii)
much better than commonly used masking techniques for input perturbation based
interpretability techniques like LIME. Thus, layer masking is able to close
the interpretability gap between CNNs and transformers, and even make CNNs
more interpretable in many cases.

##  1 Introduction

Figure 1: An outline of layer masking, our proposed method, for a
convolutional layer. The image is first masked and then padded using neighbor
padding. The convolutional layer then acts on the padded image, and a maxpool
of the same kernel size and stride acts on the mask. These are then propagated
forward through the CNN. The mask boundary is highlighted in the padded image
for illustrative purposes.

While deep learning methods have become extremely successful in solving many
computer vision tasks, they are generally opaque, and they do not admit easy
debugging of errors. Many novel interpretability methods have been developed
in recent years which attempt to analyze the reasons behind a model’s
predictions. In particular, it is natural to analyze the dependence of the
model prediction on its input by perturbing parts of its input and observing
corresponding changes in the output [18, 14, 29]. Common perturbations include
adding Gaussian noise, Gaussian blurring, replacing with a baseline color,
etc. However, many of these perturbation methods come with certain downsides.
_Partial perturbations_ , like Gaussian noise or blurring, attempt to slightly
corrupt parts of the image , while still preserving much of the information
present in those parts. While this has the advantage of not changing the input
distribution drastically, we can only measure the local sensitivity of the
model - if the model were to be robust to these perturbations, it would be
locally insensitive to perturbations on certain parts of the image but it
might still rely heavily on them for its prediction [21, 25]. _Full
perturbation_ methods remove the parts completely, and replace it with a
baseline color like black or grey. In discrete domains like natural language,
this is often the most popular method, as it is easy to remove words from the
input [15]. In vision, however, this creates a large shift in input
distribution, leading the model to perform poorly on such inputs. In
particular, it has been shown [10, 23, 24] that CNNs can be sensitive to the
mask patterns and the baseline color, focusing more on the mask than the
unmasked input.

In recent work [17, 20, 7], it has been observed that vision transformers [5]
are highly robust to many kinds of large magnitude input perturbations like
occlusions and domain shifts, maintaining upto 60% accuracy on ImageNet even
if 80% of the input is randomly blacked out. Jain et al [10] argue that this
property can make interpretability methods based on full perturbation
especially effective for transformers. They further propose to simply drop
tokens corresponding to masked out input parts instead of blacking out or
greying out image portions, just like dropping BPE tokens in a transformer-
based language model. This would make the transformer model completely
insensitive to choice of baseline color and the masking patterns.

Motivated by the same intuition, we devise a new masking technique for CNNs
which mitigates the drawbacks of full perturbation to a large extent, which we
call layer masking. Layer masking (as depicted in Fig. 1) works by running the
CNN only on the unmasked portion of the image, thus avoiding any large
distribution shift. This is done by carefully masking and padding the input of
each layer to make the model focus only on the unmasked input regions. Using
this technique, we are able to randomly remove upto 50 % of the input to a
ResNet-50 (in the form of $16\times 16$ sized patches) while maintaining the
top-1 accuracy on ImageNet over 70%, which is comparable to the robustness of
vision transformers to patch occlusion. In addition, layer masking operates at
the pixel level and is thus much more flexible than token dropping for vision
transformers which only acts on a patch level. We also find that LIME [18]
scores obtained using our masking method are more aligned with the class
object around 80% of the time as compared to simply blacking or greying out
the masked portion.

##  2 Related Work

Many interpretability methods designed for computer vision [18, 14, 29, 6]
analyze model predictions by implicitly relying on the ability to remove
features from the input. Frequently, the removed input features are either
masked out and replaced by a “neutral” baseline color like black or grey, or
perturbed partially by blurring or adding Gaussian noise. It has been shown in
[24, 23] that many of these baseline colors are not really neutral, and
interpretability methods which rely on this notion can often be quite
sensitive to choice of baselines. Secondly, the models can also be sensitive
to the shape of the masking patterns, for example, when randomly masking out
16 $\times$ 16 sized image patches in, ResNets are more likely to predict that
the image is a maze or crossword [10]. While partial perturbation methods like
adding Gaussian noise may not have similar issues, they can be misleading when
the model is insensitive to such perturbations and its output doesn’t change
significantly [21, 25]. To solve the problem of distribution shift created due
to the masking patterns, it was proposed in [9] to retrain the model with the
input perturbed according to the patterns utilized in the interpretability
methods. While this indeed solves the problem, there are two downsides: (1)
the model that is interpreted is not the original model and only a surrogate,
and thus may not be very useful for interpretability, and (2) retraining the
model with the masking pattern added is expensive. We also have methods like
inpainting the masked region using a deep generative model [3]. While this can
produce natural images and removes all the original information, it requires
training a generative model which can lead to information leakage - for
example, if a dog’s snout was masked out using this method, the generative
model may regenerate the dog’s snout again if that is the most likely
completion.

There have also been a cluster of papers [17, 20, 7, 10] showing that vision
transformers can be very robust (especially compared to CNNs) to many kinds of
perturbations including occlusions, patch permutation, adversarial
perturbations, distribution shift, etc. In connection to this, it was proposed
in [10] that interpretability methods on vision transformers are significantly
less sensitive to masking patterns and baseline colors, and also, it is
possible to simply drop the necessary patch tokens rather than substituting
with a baseline color. In line with this work, we attempt to provide a similar
masking method for CNNs which performs as well or better than token dropping
for ViTs.

In patch attacks [2], the adversary is allowed full control over a small
contiguous portion of the image. These attacks are of special interest because
they are physically realizable[11, 27] and can be a realistic way to attack
neural networks. Defenses against patch attacks [12, 19, 13, 28] use image
masking to try and mask out the portion of the image in control of the
adversary. In some of these works, the model is retrained on the images with
random features masked out to mitigate any distribution shift, as in [9],
while in others there is no retraining. We believe that our work can be
complementary to this line of work as it can improve the accuracy of the
robust classifier without needing to retrain the model.

##  3 Proposed Method: Layer Masking

###  3.1 Motivation

We design a novel feature masking technique for convolutional neural nets
(CNNs) which we call _layer masking_. We are given a model, an input image,
and a mask for the input image, and we aim to compute the model output such
that (1) it doesn’t depend on the masked out portion of the input and (2) it
only depends on the unmasked portion of the input, and not the mask pattern
itself.

Modern CNNs primarily consist of convolutional layers, along with other layers
like batch normalization, max pooling, average pooling, ReLU activations, etc.
We can categorize these layers according to the size of their receptive
fields. Layers with small receptive fields include convolutional layers, max-
pooling and average pooling layers with kernel size much smaller than the size
of the image. Fully connected layers, on the other hand, have a large
receptive fields as each output depends on all the inputs. Layers with a small
receptive field are in general more interpretable because they have fewer
parameters and are implicitly hierarchical: for example, a stack of
convolutional layers with a small kernel size processes local information
first and then progressively expands its receptive field to encompass the
whole image. We exploit this structure by devising an algorithm which
carefully masks the input and output at each layer with small receptive field
such that information loss and artifacts created by the masking procedure at
each step is minimal. We propagate both the input and the mask at each layer
so as to simulate running a CNN on an irregularly shaped input corresponding
to the unmasked input features, rather than substituting the masked inputs
with a baseline color. We are careful, however, to not propagate forward any
information in the masked out input regions.

Let the input to a convolutional layer with small receptive field be
${\bm{x}}\in{\mathbb{R}}^{c\times n\times n}$ with output
${\bm{y}}\in{\mathbb{R}}^{c^{\prime}\times n^{\prime}\times n^{\prime}}$ and
binary input mask ${\bm{m}}\in\\{0,1\\}^{n\times n}$ (${\bm{m}}[u,v]=1$
implies that cell $(u,v)$ is unmasked, else it is masked out). Each element of
the output of this layer with kernel size $k\times k$ depends on at most
$k^{2}$ input values. These input values may either be all masked, all
unmasked or partially masked and unmasked (when the convolution is over the
mask edge), depending on the values of ${\bm{m}}$ over the receptive field.

It is clear that our masking procedure should propagate forward the outputs
which only depend on the unmasked input, and discard those outputs which
depend only on the masked portion. However, it is not immediately obvious how
to handle the outputs from the convolutions over the mask edge. The challenge
here is that edge convolutions contain valuable information about the edges,
and if we discard them at each layer, the unmasked portion of the image can
quickly vanish to zero. Thus, we choose to propagate forward the edge
convolutions. However, there is the danger of them distorting the natural
distribution of the layer activations, as the output unavoidably depends on
the masked out region which is filled with zeros. For example, in the third
figure (bottom row) of Fig. 2, we see a slice of the activations obtained
after applying the 1st residual block of ResNet-50 on the image with the
central square region masked out at every layer, but including all the edge
convolutions in the output. We see that the convolutions at the top edge of
the mask result in a brighter top edge which indicates high activations. This
is undesirable since this is an artifact created due to the masking method. We
hypothesize that this is because the abrupt transition between the unmasked
input and zeros trigger the filters sensitive to edges, thus creating a large
activation.

To mitigate this issue, rather than just fill the masked out portion with
zeros, we pad the unmasked portion using a variant of replication padding we
call neighbor padding. Specifically, we iteratively assign the masked input
cells adjacent to the mask edge with the average value of its immediate non-
zero neighbors. This process is continued till the width of the padding is at
least $k$, the kernel size of the layer. In Fig. 3, we see that after the
cells near the edge are progressively filled using the values of its
neighbors, the resultant image looks very natural and there is no sharp
discontinuity near the edge. We find experimentally that this works much
better than padding with zeros (see Tab. 3).

Algorithm 1 Neighbor padding algorithm

Input: Input to be padded ${\bm{x}}$, Mask ${\bm{m}}$, padding width $k$

Output: Padded input ${\bm{x}}^{\prime}$

Initialize ${\bm{x}}^{\prime}\leftarrow{\bm{x}}\odot{\bm{m}}$,
$\epsilon\leftarrow 10^{-8}$

Initialize ${\bm{f}}\leftarrow{\bm{1}}_{3\times 3}$, a $3\times 3$ filter
filled with ones

for $i=1$ to $k$ do

${\bm{n}}\leftarrow{\bm{x}}^{\prime}*{\bm{f}}$ // Numerator of the neighbor
average

${\bm{d}}\leftarrow{\bm{m}}*{\bm{f}}$ // Denominator of the neighbor average

${\bm{e}}\leftarrow(1-{\bm{m}})\odot{\bm{n}}/({\bm{d}}+\epsilon)$ // Fill
masked inputs

${\bm{x}}^{\prime}\leftarrow{\bm{x}}^{\prime}+{\bm{e}}$

${\bm{m}}\leftarrow\min({\bm{1}},{\bm{m}}+{\bm{d}})$ // Update masks

end for

We also have to propagate the masks forward, such that for the output of any
layer, the propagated mask at that layer is of the same shape as the output
and indicates which output values need to be masked out by the following
layers. Since edge convolutions are not discarded at any step, the propagated
mask must contain 1 for all output cells which depends on the unmasked portion
of the input, and 0 everywhere else.

###  3.2 Formal Description

We now describe our method more formally. Suppose we are given a CNN $f$ which
is structured like a directed acyclic graph. Each node of the DAG represents a
layer or operation which acts on the outputs of the nodes with which it has an
incoming edge. We replace each layer with its masking version (subscripted
with $m$) which acts on an input-mask pair. Let $g_{k}$ be a layer with
receptive field of size $k$. Then, we define its masking version:

$g_{k,m}({\bm{x}},{\bm{m}})=(g_{k}(\text{Pad}_{k}({\bm{x}},{\bm{m}})),\text{MaxPool}_{k}({\bm{m}})))$

In this equation, $g_{k}$ could be any convolutional or pooling layer with
kernel size $k$ and some stride $s$. $\text{MaxPool}_{k}$ is a max pooling
layer with the same kernel size and stride as $g_{k}$.
$\text{Pad}_{k}({\bm{x}},{\bm{m}})$ is a function which neighbor pads
${\bm{x}}\odot{\bm{m}}$ with padding width $k$ (described in Algorithm 1).
Here, $\odot$ is the Hadamard product (with suitable broadcasting), and $*$ is
convolution with zero padding. The max pool layer ensures that the output
masks contains a 1 for all convolutions where even a single input was
unmasked.

Layers which act independently on each element (like ReLU and BatchNorm) can
be considered to be a special case of the above with $k=1$. In this case, the
above equation is greatly simplified and becomes:

$g_{m}({\bm{x}},{\bm{m}})=(g({\bm{x}}\odot{\bm{m}}),{\bm{m}})$

In models like ResNet which use residual connections, two sets of activations
can be added together. We can add two input - mask pairs as:

$({\bm{x}}_{1},{\bm{m}}_{1})+({\bm{x}}_{2},{\bm{m}}_{2})=(({\bm{x}}_{1}+{\bm{x}}_{2})\odot({\bm{m}}_{1}\odot{\bm{m}}_{2}),{\bm{m}}_{1}\odot{\bm{m}}_{2})$

Note that we lose some information here by taking the Hadamard product of the
masks, but this is negligible in practice. For fully connected layers, there
are hardly any output values which depend on only the masked or unmasked
portion of the image. Thus, we allow the fully connected layer to operate
normally on its input, replacing the masked portions of the input with a zero
baseline value. If $h$ is a fully connected layer, then we define its masking
version:

$h_{m}({\bm{x}},{\bm{m}})=h({\bm{x}}\odot{\bm{m}})$

Any layers after a fully connected layer act on the input as normal. We can
now create a new model $f_{m}$ which has the same DAG structure of the
original model $f$, except each layer $g^{i}$ or $h^{i}$ has been replaced
with the corresponding masking version $g^{i}_{m}$ or $h^{i}_{m}$. $f_{m}$
acts on an image - mask pair and produces an output which depends only on the
unmasked portion of the image.

Figure 2: Top row: A dog image, sample activations after the first residual
block, and mask to be applied to the image. Bottom row: The same activations
when using grey-out masking, layer masking with neighbor padding (our method),
and layer masking without neighbor padding (using zero padding). Neighbor
padding helps in eliminating undesirable edge artifacts encountered in zero
padding and greying out. Layer masking completely zeros out the masked out
region unlike greying out which has non-zero values after a few layers Figure
3: A visual depiction of neighbor padding on a part of the image as $k$
increases. The grey line is the mask edge (added in for illustrative
purposes), the cells near the edge are filled progressively with the average
of their neighbors’ values

##  4 Experiments

We perform our experiments on CNNs (ResNet-50 [8]) and vision transformers
(ViT-B/16 [5]) pretrained on ImageNet-1000 [4]. The datasets we use are thus
also subsets of ImageNet. We analyze and compare the following input masking
methods: (1) Layer masking, as described in Sec. 3 (only for CNNs), (2) Token
dropping, as described in [10] (only for ViTs), (3) Grey-out: Replace the
masked out input with a grey baseline color equal to the ImageNet mean and (4)
Black-out: Replace the masked out input with a black baseline color.

Figure 4: Segmenting (a) an image of a dog using (b) 16 $\times$ 16 patches,
(c) SLIC and (d) Quickshift

We also use the following segmentation algorithms to extract the segments from
the image (see Fig. 4 for a pictorial representation of these) (1) Square
patches: Segment the $224\times 224$ image into smaller square patches
(typically $16\times 16$), (2) SLIC superpixels [1]: Segment the image by
clustering pixels based on similarity in the color space and local proximity,
and (3) Quickshift [26]: Segment the image by first computing a kernel density
estimate of the pixel distribution, and then connect each pixel to its nearest
neighbor with higher density. We tune hyperparameters for these algorithms
such that they divide the image into approximately the same number of
segments.

(a) Accuracy

(b) Class entropy

Figure 5: Changes in model accuracy and class entropy vs fraction of 16
$\times$ 16 patches(1st row), Quickshift segments (2nd row) and SLIC
segments(3rd row) masked out in (1st column) random order, (2nd column) most
salient first, and (3rd column) least salient first. ResNets with layer
masking are clearly less sensitive to segment removal, especially removal of
$16\times 16$ patches. They are also unbiased by the mask shape, unlike
ResNets with greyout or blackout masking which predict classes like maze and
crossword when $16\times 16$ patches are blacked out or greyed out. ResNets
with layer masking match or exceed robustness of ViTs to patch occlusions.
Note: Quickshift produces many very small segments which are not very salient,
thus when removing the half of the most salient quickshift segments, most of
the image is removed, which is why accuracy and class entropy falls more
quickly for quickshift.

###  4.1 Segment ablation experiments

To quantify the effect of feature masking methods on the model predictions, we
study the behavior of the model when varying parts of the input are masked
out. We first segment the image using a segmentation algorithm. Then, we
analyze how the model output changes when more and more segments are masked
out using the above masking techniques. We characterize the model behavior
using 4 metrics: accuracy, class entropy (defined as entropy of
$p_{f}(y)=\mathbb{E}_{{\bm{x}}\in D}[1[f({\bm{x}})=y]]$), WordNet similarity
[16] between predictions and true labels, and fraction of unchanged
predictions. The WordNet similarity measures how similar the model predictions
are to the true labels on a scale from 0 to 1 in place of a binary hit / miss.
The fraction of unchanged predictions is a measure of the number of
predictions changed by masking out parts of the input image.

In line with previous work [17, 10], we mask out these segments in three ways:
(1) Randomly, (2) Most salient first, (3) Least salient first. We then
evaluate the metrics listed above, and plot them as a function of the fraction
of segments masked out. The accuracy and class entropy plots are shown in Fig.
5, while the rest can be found in the appendix. For the sake of clarity, we
report only the best performing baseline for vision transformers, which is
token dropping [10] when the segments are $16\times 16$ patches, but greyout
masking for other segmentation methods.

Computing the saliency scores: We select 5000 random images and saliency maps
from Salient ImageNet [22] for our experiment. Each saliency map
$m\in{\mathbb{R}}^{d\times d}$ is a pixel level saliency attribution array
where $0\leq m[i,j]\leq 1$ denotes the importance of the $(i,j)$th pixel to
predicting the ImageNet class - the higher the value, the more important the
pixel is. Each image is first segmented using a given segmentation algorithm.
We then compute saliency scores for each segment by summing up the saliency
values for all pixels in that segment.

We expect that an ideal masking technique should affect the model output
minimally, under the constraint that the model does not rely on information in
the parts of the image which are masked out. It should also not be biased by
the masking pattern to predict certain classes more than others. We measure
the model prediction quality using accuracy, WordNet similarity and change in
model predictions, and the extent of masking bias using class entropy - the
higher the class entropy, the less biased the masking.

We observe that, consistent with [17, 10], vision transformers are more robust
to image occlusion as compared to ResNets when greying out or blacking out the
image. Specifically, we find that transformers are much more robust to
perturbations which are aligned with the patch boundaries of the vision
transformer (like occluding $16\times 16$ sized patches for ViT-B/16), but are
more sensitive to occlusions of irregularly shaped regions, as found in SLIC
or quickshift segments. However, if layer masking is used, ResNets can be much
more robust to such occlusions irrespective of its shape, exceeding robustness
of transformers in many cases. The class entropy plots show that layer masking
does not bias the model towards predicting certain classes. Greying or
blacking out segments however, makes the model predict certain classes more
than others - for example, greying out $16\times 16$ patches causes the model
to predict ImageNet classes like ‘maze’ or ‘crossword’ which resemble the
masking pattern [10]. We thus find that layer masking is superior to greying
out or blacking out across all cases.

###  4.2 Impact on interpretability methods

Many interpretability methods [18, 14, 29] rely on the notion of a “neutral”
masking technique in order to evaluate the influence of different parts of the
image on the model prediction. This is generally not true and the masking
technique can have a significant impact on the importance scores given by the
interpretability method. We investigate the effect of masking methods on
interpretability in this section using the example of LIME.

####  4.2.1 Visual inspection

(a) A cat and a mouse

(b) A dog and a cat

Figure 6: Visualization of LIME scores for the top two predictions of
ResNet-50 on two sample images. Columns correspond to the masking techniques
(blacking out, greying out, and layer masking), rows are the top 2
predictions. For (a), the top two predictions are American black bear and
Groenendael. For (b), they are Bernese mountain dog and tabby cat. Green
regions contribute to the prediction, red regions detract from the prediction.

We first visually illustrate the impact of masking technique on LIME using a
couple of examples in Fig. 6. For the LIME explanations, the top 10 segments
with the highest magnitude LIME score are highlighted. Red segments detract
from the prediction, and a green segments contribute to the prediction. In the
image of a cat and a mouse, the model misclassifies the black cat as an
‘American black bear’ or a ‘groenendael’ (a black dog breed). But LIME assigns
a highly negative score (in red) to certain regions of the black cat for the
‘groenendael’ and ‘black bear’ prediction when masking using blacking or
greying out, which seems incorrect. Similarly, in the image of a cat and a
dog, LIME with greyout or blackout masking assigns parts of the cat high
positive scores (in green) and parts of the dog negative scores (in red) in
the explanation of the prediction of Bernese mountain dog, and vice versa for
the explanation of the tabby prediction. The third column containing the LIME
explanations obtained using layer masking is much more aligned with human
intuition.

####  4.2.2 Quantitative evaluation

| Blackout | Greyout | Layer masking
---|---|---|---
ResNet-50 |  |  |
Quickshift | 0.1153 | 0.2211 | 0.3042
SLIC | 0.1235 | 0.2186 | 0.3120
16 x 16 patches | 0.0119 | 0.0675 | 0.2253
ViT/B-16 |  |  | Token dropping
Quickshift | 0.1693 | 0.2417 | 0.0444
SLIC | 0.1748 | 0.2446 | 0.0539
16 x 16 patches | 0.1633 | 0.1678 | 0.1430
Table 1: Mean alignment scores for LIME computed over 500 randomly sampled images from Pixel ImageNet [30] | Blackout | Greyout | Layer masking
---|---|---|---
ResNet-50 |  |  |
Quickshift | 0.106 | 0.228 | 0.666
SLIC | 0.078 | 0.206 | 0.716
16 x 16 patches | 0.040 | 0.146 | 0.814
ViT/B-16 |  |  | Token dropping
Quickshift | 0.292 | 0.594 | 0.114
SLIC | 0.270 | 0.614 | 0.116
16 x 16 patches | 0.366 | 0.396 | 0.238
Table 2: Mean win rate for LIME computed over 500 randomly sampled images from
Pixel ImageNet

To get a better sense of the effect of masking techniques on LIME, we perform
a quantitative evaluation. For each prediction, LIME assigns an importance
score $s_{i}$ to each image segment $i$, where a positive score of high
magnitude means that the segment contributes highly to the prediction, but a
negative score means that the segment detracts from the prediction. To
evaluate explanations quantitatively in the absence of a reliable ground
truth, we use the images and segmentation masks from Pixel ImageNet [30] as a
substitute for the “correct” explanation. Given a segmentation mask
${\bm{m}}\in[0,1]^{d\times d}$ for an image of dimension $d$, we derive the
“ground truth” ${\bm{g}}$ for the explanation such that
${\bm{g}}_{i}=\sum_{(u,v)\in\text{patch }i}\left({\bm{m}}[u,v]-m_{avg}\right)$
where $m_{avg}$ is the mean of the segmentation mask. We can then measure how
aligned the explanations are with the ground truth by computing the _alignment
score_ , which is the cosine similarity between $g_{i}$ and $s_{i}$, or
$\cos({\bm{g}},{\bm{s}})=\frac{\sum_{i}g_{i}s_{i}}{\sqrt{(\sum_{i}g_{i}^{2})(\sum_{i}s_{i}^{2})}}$.
The alignment score will be $1$ if the LIME explanation ${\bm{s}}$ is
perfectly aligned with ${\bm{g}}$, and $-1$ if it is completely misaligned.

We measure alignment scores for different masking methods for ResNet-50 and
ViT over 500 images from Pixel ImageNet and report them in Tab. 1. In these
500 images, the correct class is within the top 3 predictions of the model. In
general, ViTs have better scores for greyout and blackout than the
corresponding scores for ResNet-50, especially for $16\times 16$ sized
patches. We find that token dropping in vision transformers is not the best
solution - although it performs somewhat well if the segments are $16\times
16$ patches, it fails when the segments are irregularly shaped. This is
because token dropping operates at the patch level and has to discard more
patches and more information than greying out or blacking out in order for it
to not rely on any masked out information. However, layer masking does not
suffer the same problems. We find that layer masking is much superior to
blacking or greying out the input across different segmentation algorithms. In
fact, the alignment scores for ResNet-50 with layer masking are better than
the best numbers for ViTs, which suggests that ResNets can be more
interpretable than ViTs when combined with layer masking.

However, note that the alignment scores can vary quite a lot over many images.
For example, there might be some images on which LIME works well, and others
where it fails, leading to a large variation in alignment scores. Therefore,
it might be misleading to just report the average alignment score. We also
calculate the _win rate_ for each masking technique, which is the fraction of
images for which that particular masking technique had the best alignment
score compared to the other techniques in Tab. 2. The win rates sum up to 1 in
each row. We find that trends are broadly similar to those in Tab. 1. Layer
masking wins over regular black or grey masking 65% to 80% of the time.

###  4.3 Ablation study

To further investigate the effect of layer masking and neighbor padding on
model behavior, we construct 3 variants of layer masking: (a) With zero
padding instead of neighbor padding (b) Masking and padding only the first two
residual blocks, (c) Masking and padding only the first convolutional layer,
ReLU and BatchNorm layer

Using a similar setup as in Sec. 4.1, we we compute the area under curve (AUC)
for each plot of metric vs fraction of segments dropped. The AUC values are
averaged over different segmentation algorithms (SLIC, quickshift, etc) and
masking orders (random, salient first, etc)(refer Tab. 3).

We find that both neighbor padding and masking all layers are important to the
masking technique. Layer masking with zero padding is still better than
blackout or greyout, but much worse than with neighbor padding. Layer masking
only the first two residual blocks is also inferior to masking through all
layers, but we find that there are diminishing returns, as we are able to
obtain much of the improvement by masking only half of the layers.

|  Accuracy |  Class |  Wordnet |  Unchanged
---|---|---|---|---
|  |  Entropy |  Similarity |  Predictions
Blackout |  0.3881 |  4.6473 |  0.6930 |  0.4094
Greyout |  0.4398 |  4.9408 |  0.7167 |  0.4636
Layer masking: |  |  |  |
On all layers |  0.5604 |  5.6021 |  0.7881 |  0.5907
On 1st and 2nd residual blocks |  0.5103 |  5.0962 |  0.7616 |  0.5391
Zero padding |  0.4502 |  5.0388 |  0.7262 |  0.4747

Table 3: Average AUC for different variants of layer masking alongside the
black out and grey out baselines (model: ResNet-50). Higher the better

###  4.4 Other interesting properties of layer masking

In this section, we identify some more properties of layer masking that are
important for model interpretability.

Linearity in masking: Consider a model equipped with a masking technique
$f_{m}$ which acts on an input - mask pair $({\bm{x}},{\bm{m}})$ and returns
an output ${\bm{y}}$ which depends only on the unmasked parts of the input.
Then, we say that the model $f_{m}$ is linear in masking if
$f_{m}({\bm{x}},{\bm{m}}_{1}+{\bm{m}}_{2})=f_{m}({\bm{x}},{\bm{m}}_{1})+f_{m}({\bm{x}},{\bm{m}}_{2})$
for any two binary masks ${\bm{m}}_{1},{\bm{m}}_{2}$ such that
${\bm{m}}_{1}\cdot{\bm{m}}_{2}=0$. This property is useful for
interpretability methods like LIME which train a linear model on
$({\bm{m}},{\bm{y}})$ pairs and use its weights to explain the model
prediction. Modern vision models like CNNs and Vision Transformers are non-
linear and include cross-interactions between features in ${\bm{m}}_{1}$ and
${\bm{m}}_{2}$. Thus, it is not possible to design a perfectly linear masking
technique for these architectures, which means that only approximate linearity
is possible. However, we can attempt to design more linear masking methods for
each model architecture, and thus obtain more interpretable masking
techniques.

We measure linearity by sampling random images from ImageNet and dividing it
into $N$ smaller square patches. We can then compute the cosine similarity
between $f({\bm{x}})$ and $\sum_{i=1}^{N}f_{m}({\bm{x}},{\bm{m}}_{i})$ where
${\bm{m}}_{i}$ corresponds to patch $i$ (Tab. 4). We find that layer masking
is much more linear as compared to greying out or blacking out pixels, and in
general, ResNet masking methods are more linear than corresponding methods for
ViTs. Because the attention heads in ViTs introduce a lot of cross terms right
from the beginning, including cross terms between distant patches, linearity
in vision transformer masking is much lower than CNN masking.

We also find that in layer masking,
$\mathbb{E}_{\bm{x}}\|f_{m}({\bm{x}},{\bm{m}})\|$ scales linearly with
$|{\bm{m}}|$. We test this by measuring the magnitude of
$f_{m}({\bm{x}},{\bm{m}})$ with ${\bm{m}}$ as a mask for square patches of
side length $n$, so that $\|{\bm{m}}\|\propto n^{2}$. We observe in Fig. 8
that layer masking closely tracks the $n^{2}$ curve, which implies that
$\mathbb{E}_{\bm{x}}\|f_{m}({\bm{x}},{\bm{m}})\|$ scales almost linearly with
$\|{\bm{m}}\|$ for layer masking. However, the magnitude for ViT features
remain approximately constant.

Avoidance of output collapse: As the fraction of masked input approaches 1, it
is desirable to avoid the model output collapsing to the same vector and thus
not being sensitive enough to the unmasked features. To test for this, we take
two random images ${\bm{x}}_{1}$ and ${\bm{x}}_{2}$ of size $224\times 224$
and compute the cosine similarity between their model features,
$c=\cos(f({\bm{x}}_{1}),f({\bm{x}}_{2}))$. Then, these images are resized to a
smaller size $n$, and padded with zeros to recover the original size . We now
have images ${\bm{x}}_{1,n}$ and ${\bm{x}}_{2,n}$ of size $224\times 224$ and
a mask of the same shape ${\bm{m}}_{n}$ which is 1 for a region of size
$n\times n$ and 0 elsewhere. We then measure the cosine similarity between
${\bm{x}}_{1,n}$ and ${\bm{x}}_{2,n}$ as
$c_{n}=\cos(f_{m}({\bm{x}}_{1,n},{\bm{m}}_{n}),f_{m}({\bm{x}}_{2,n},{\bm{m}}_{n}))$
and plot $\mathbb{E}_{{\bm{x}}_{1},{\bm{x}}_{2}}[c_{n}-c]$ as function of $n$
in Fig. 7. We clearly see that as the image size is decreased, the cosine
similarity changes much more for greyout or blackout as compared to layer
masking for ResNet-50 or token dropping for ViTs.

Figure 7: Average difference in cosine similarity vs image size. Since model
features of ViTs can be negative unlike ResNet-50, cosine similarity can vary
from -1 to 1  Figure 8: Mean magnitude of output feature vectors vs image size

| ResNet-50 | ViT-B/16
---|---|---
Patch size |  Blackout |  Greyout |  Layer masking |  Blackout |  Greyout |  Token dropping
112 |  0.7975 |  0.8284 |  0.9485 |  0.6707 |  0.7063 |  0.7043
56 |  0.5124 |  0.5842 |  0.8310 |  0.2202 |  0.2506 |  0.1929
32 |  0.4282 |  0.4878 |  0.7094 |  0.1377 |  0.1365 |  0.1426
16 |  0.3848 |  0.4371 |  0.6490 |  0.0912 |  0.0876 |  0.0877

Table 4: Average cosine similarity between image features and their linear
approximation

##  5 Conclusion

In this paper, we have presented a new masking technique such that the model
output is both (a) perfectly insensitive to the masked out portion of the
input and (b) only focused on the unmasked input and not the masking pattern.
We find that layer masking can make CNNs like ResNets very robust to input
removal without retraining - comparable to vision transformers - especially
when the masking patterns can get confused with output classes like patch
occlusions. We further find that LIME scores obtained using layer masking are
better than LIME scores with blacking or greying out for up to 80% of the
images in Pixel ImageNet. We propose that this technique can be of great use
in interpretability techniques like LIME which rely on the ability to remove
features from the input without any major distribution shift.

## References

  * [1] Radhakrishna Achanta, Appu Shaji, Kevin Smith, Aurelien Lucchi, Pascal Fua, and Sabine Süsstrunk.  Slic superpixels compared to state-of-the-art superpixel methods.  IEEE Transactions on Pattern Analysis and Machine Intelligence, 34(11):2274–2282, 2012.
  * [2] Tom B Brown, Dandelion Mané, Aurko Roy, Martín Abadi, and Justin Gilmer.  Adversarial patch.  arXiv preprint arXiv:1712.09665, 2017.
  * [3] Chun-Hao Chang, Elliot Creager, Anna Goldenberg, and David Duvenaud.  Explaining image classifiers by counterfactual generation.  In International Conference on Learning Representations, 2019.
  * [4] Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei.  Imagenet: A large-scale hierarchical image database.  In 2009 IEEE conference on computer vision and pattern recognition, pages 248–255. Ieee, 2009.
  * [5] Alexey Dosovitskiy, Lucas Beyer, Alexander Kolesnikov, Dirk Weissenborn, Xiaohua Zhai, Thomas Unterthiner, Mostafa Dehghani, Matthias Minderer, Georg Heigold, Sylvain Gelly, Jakob Uszkoreit, and Neil Houlsby.  An image is worth 16x16 words: Transformers for image recognition at scale.  In International Conference on Learning Representations, 2021.
  * [6] Ruth C Fong and Andrea Vedaldi.  Interpretable explanations of black boxes by meaningful perturbation.  In Proceedings of the IEEE international conference on computer vision, pages 3429–3437, 2017.
  * [7] Jindong Gu, Volker Tresp, and Yao Qin.  Are vision transformers robust to patch perturbations?  In Shai Avidan, Gabriel Brostow, Moustapha Cissé, Giovanni Maria Farinella, and Tal Hassner, editors, Computer Vision – ECCV 2022, pages 404–421, Cham, 2022. Springer Nature Switzerland.
  * [8] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  CoRR, abs/1512.03385, 2015.
  * [9] Sara Hooker, Dumitru Erhan, Pieter-Jan Kindermans, and Been Kim.  A benchmark for interpretability methods in deep neural networks.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019.
  * [10] Saachi Jain, Hadi Salman, Eric Wong, Pengchuan Zhang, Vibhav Vineet, Sai Vemprala, and Aleksander Madry.  Missingness bias in model debugging.  In International Conference on Learning Representations, 2022.
  * [11] Mark Lee and Zico Kolter.  On physical adversarial patches for object detection.  arXiv preprint arXiv:1906.11897, 2019.
  * [12] Alexander Levine and Soheil Feizi.  (de)randomized smoothing for certifiable defense against patch attacks.  In Proceedings of the 34th International Conference on Neural Information Processing Systems, NIPS’20, Red Hook, NY, USA, 2020. Curran Associates Inc.
  * [13] Jiang Liu, Alexander Levine, Chun Pong Lau, Rama Chellappa, and Soheil Feizi.  Segment and complete: Defending object detectors against adversarial patch attacks with robust patch detection.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 14973–14982, 2022.
  * [14] Scott M Lundberg and Su-In Lee.  A unified approach to interpreting model predictions.  In I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett, editors, Advances in Neural Information Processing Systems 30, pages 4765–4774. Curran Associates, Inc., 2017.
  * [15] Dina Mardaoui and Damien Garreau.  An analysis of lime for text data.  In Arindam Banerjee and Kenji Fukumizu, editors, Proceedings of The 24th International Conference on Artificial Intelligence and Statistics, volume 130 of Proceedings of Machine Learning Research, pages 3493–3501. PMLR, 13–15 Apr 2021.
  * [16] George A. Miller.  WordNet: A lexical database for English.  In Human Language Technology: Proceedings of a Workshop held at Plainsboro, New Jersey, March 8-11, 1994, 1994.
  * [17] Muzammal Naseer, Kanchana Ranasinghe, Salman Khan, Munawar Hayat, Fahad Khan, and Ming-Hsuan Yang.  Intriguing properties of vision transformers.  In A. Beygelzimer, Y. Dauphin, P. Liang, and J. Wortman Vaughan, editors, Advances in Neural Information Processing Systems, 2021.
  * [18] Marco Tulio Ribeiro, Sameer Singh, and Carlos Guestrin.  ”why should I trust you?”: Explaining the predictions of any classifier.  In Proceedings of the 22nd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, San Francisco, CA, USA, August 13-17, 2016, pages 1135–1144, 2016.
  * [19] Hadi Salman, Saachi Jain, Eric Wong, and Aleksander Madry.  Certified patch robustness via smoothed vision transformers.  CoRR, abs/2110.07719, 2021.
  * [20] Rulin Shao, Zhouxing Shi, Jinfeng Yi, Pin-Yu Chen, and Cho-Jui Hsieh.  On the adversarial robustness of vision transformers, 2021.
  * [21] Avanti Shrikumar, Peyton Greenside, and Anshul Kundaje.  Learning important features through propagating activation differences.  In Doina Precup and Yee Whye Teh, editors, Proceedings of the 34th International Conference on Machine Learning, volume 70 of Proceedings of Machine Learning Research, pages 3145–3153. PMLR, 06–11 Aug 2017\.
  * [22] Sahil Singla and Soheil Feizi.  Salient imagenet: How to discover spurious features in deep learning?  In International Conference on Learning Representations, 2022.
  * [23] Pascal Sturmfels, Scott Lundberg, and Su-In Lee.  Visualizing the impact of feature attribution baselines.  Distill, 2020.  https://distill.pub/2020/attribution-baselines.
  * [24] Mukund Sundararajan and Ankur Taly.  A note about: Local explanation methods for deep neural networks lack sensitivity to parameter values.  CoRR, abs/1806.04205, 2018.
  * [25] Mukund Sundararajan, Ankur Taly, and Qiqi Yan.  Axiomatic attribution for deep networks.  In International conference on machine learning, pages 3319–3328. PMLR, 2017.
  * [26] Andrea Vedaldi and Stefano Soatto.  Quick shift and kernel methods for mode seeking.  In David Forsyth, Philip Torr, and Andrew Zisserman, editors, Computer Vision – ECCV 2008, pages 705–718, Berlin, Heidelberg, 2008. Springer Berlin Heidelberg.
  * [27] Zuxuan Wu, Ser-Nam Lim, Larry S Davis, and Tom Goldstein.  Making an invisibility cloak: Real world adversarial attacks on object detectors.  In European Conference on Computer Vision, pages 1–17. Springer, 2020.
  * [28] Maksym Yatsura, Kaspar Sakmann, N Grace Hua, Matthias Hein, and Jan Hendrik Metzen.  Certified defences against adversarial patch attacks on semantic segmentation.  arXiv preprint arXiv:2209.05980, 2022.
  * [29] Matthew D. Zeiler and Rob Fergus.  Visualizing and understanding convolutional networks.  CoRR, abs/1311.2901, 2013.
  * [30] Shiyin Zhang, Jun Hao Liew, Yunchao Wei, Shikui Wei, and Yao Zhao.  Interactive object segmentation with inside-outside guidance.  In Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition, pages 12234–12244, 2020.


# Electrostatics in a crooked nanochannel in a newly developed curvilinear
coordinate system

Xi Chen, Ke Xiao, Rui Ma, Xuezheng Cao, and Chen-Xu Wu<EMAIL_ADDRESS>Department of Physics, School of Physical Science and Technology, Xiamen
University, Xiamen 361005, People’s Republic of China

###### Abstract

Both biological and artificial nanochannels in crooked shape exhibit unusual
transportational characteristics, bringing about a challenge to the
traditional theoretical analysis of nanofluidics, partly due to their
complicated boundary description. In this paper, by developing a curvilinear
coordinate system for crooked nanochannels, we successfully solve the
electrostatic Poisson-Boltzmann equation analytically for a two-dimensional
nanochannel, with its effectiveness confirmed through numerical calculation.
The influences of the geometric profile of the nanochannel on the distribution
of electric potential, ionic concentration, and surface charge on channel
walls can be quantitatively evaluated in a facilitated way in terms of these
curvilinear coordinates. Such a technique can be widely applied to many
nanofluidic systems.

##  I Introduction

Nanochannels containing ionic solution in both biological and artificial
systems [1, 2, 3] have exhibited many specific features that are absent in
macroscopic configurations. In a nanoscale confinement, the interaction
between ions and the channel becomes comparatively significant since the
suface-to-volume ratio is high. In some biological nanochannels, their wall
usually is composed of some functional groups, which will be ionized once in
contact with electrolyte solution. Such kind of ionization occurring in the
solid-liquid interface makes the channel wall deviate from electroneutrality
and greatly influence the ionic distribution in the channel, which, as a
result, impacts electrical properties of the whole channel. Many
investigations have increasingly been focused on the effect of interfaces
between nanochannel walls and electrolyte solution [4, 5, 6].

In nanoscale, the electrostatic distribution of charged particles near a
solid-liquid interface is described as electric double layer (EDL) where
charged channel wall attracts counterions to maintain electroneutrality.
According to Gouy-Chapman theory [7, 8], ions in EDL can be regarded as point
charges dispersed in uniform dielectric solvent with constant permitivity, so
the continuum theory is still valid to describe such nanosystem. In
theoretical practice, Poisson-Boltzmann (PB) equation is applied to construct
a mean-field model for EDL and its analytical solutions can be expressed in
terms of Cartesian coordinates [9], cylindrical coordinates [10], or spherical
coordinates [11, 12], depending on the particular profile of the boundaries.

Previous theoretical investigations on nanochannels were mainly focused on
configurations with straight axis, such as cylinder [13, 14, 15], cone [16],
or two-dimensional (2D) multilayered-system [17] and funnel [18], whose axial
symmetry can be taken advantage of to simplify governing equations and their
corresponding constraint equations in a 2D or even 1D fashion. However, in
most cases when such an axial symmetry no longer exists, it becomes extremely
hard to analytically solve their governing equations. In Constantin and Siwy’s
investigation on a conical nanochannel [16], for example, electric potential
and ionic concentration were averaged over the transverse coordinates and
formulated as functions of axial coordinate. For comparison, they also carried
out a numerical simulation for the whole ionic current through the nanochannel
based on the Poisson-Nernst-Planck (PNP) equations without any assumption in
advance. When the nanochannel is in an aformensioned regular shape, the
analytical solution of PB equation is still attainable via approximation, but
in cases of broken axial symmetry, usually one has to resort to numerical
simulation instead due to the lack of analytical mathematical tool.

Recent experiments and modelings have discovered some novel phenomena like
transportational properties affected by the profile of nanochannels, e.g.,
curvature plays a significant role in the intercellular exchange between
organelles [19, 20] and the electrokinetic flow in artificial nanochannels
[21]. Unlike straight nanochannels, the distribution of charged particles in a
curved channel are governed by curvature-dependent electric potential, which
provides a method to dynamically control the directional ion transport inside,
i.e., current rectification [22, 23]. Meanwhile, in recent nanofluidic
experiments, ionic rectification was observed in a system composed of
nanosheets fabricated on flexible materials[24], which can be treated as a 2D
nanochannel with asymmetric geometric configuration, i.e. a 2D curved
nanochannel. Such kind of 2D nanochannel can also be found in nanosheets
embedded in tunable membranes, a system also showing dynamic shape effect[25,
26].

Although theoretical analyses of curvature effect have been done in
rectangular microchannels with uniform curvature [27, 28], curved nanochannels
with varying curvatures along their contours are rarely discussed. Due to the
irregularity of its geometric configuration, it is difficult to express the
boundary conditions for an electrostatic behavior like EDL based on any
traditional coordinate systems, let alone find its analytical solution. Thus
it is necessary to develop a suitable coordinate system so as to describe the
curved space inside the nanochannel in an easy manner.

In this paper, we develop a curvilinear coordinate system based on Frenet
frame, which is found to greatly simplify the representation of the space
inside a curved channel and its boundaries. Due to the difficulties in fully
solving 3D electrostatic equations analytically, we derive the analytical
solution to the PB equation governing the electrostatics inside a 2D curved
nanochannel, like the systems reported in Refs. [24, 25, 26]. The effect of
channel profile on distribution of electric potential, ionic concentration and
surface charge under different boundary conditions is also analyzed. The
numerical simulations carried out under the same conditions using finite
element method support the effectiveness of our theoretical analysis under
such a new curvilinear coordinate frame.

##  II Theoretical analysis of PB equation based on the new curvilinear
coordinate system

###  II.1 Curvilinear coordinate system

In differential geometry theory, Frenet frame is a right-handed orthogonal
unit frame moving along a canonical curve with non-zero curvature, which can
be defined by Frenet-Serret formula [29]. Numerous studies have shown its
efficacy and convenience in representing vector fields of physical properties
inside a curved space, like the electromagnetic field of plasma. For example,
for the purpose of analyzing magnetohydrodynamics of the magnetic flux in a
toroidal space, like tokamaks [30, 31], Frenet frame was utilized to establish
an orthogonal curvilinear coordinate frame for the equation of Lorentz force,
which depends on the curvature of the tube. Using Frenet frame, the central
axis of a planarly crooked nanochannel we consider in $XY$ plane can be
expressed, in terms of arc length parameter $s$, as

$\bm{r_{0}}(s)=(X(s),Y(s),0).$ (1)

The unit tangent vector $\alpha(s)$ at a point $s$ of the axial curve, as
shown in Fig. LABEL:fig1(a), is defined as

$\bm{\alpha}(s)=(X^{\prime}(s),Y^{\prime}(s),0),$ (2)

where the primes refer to the derivatives with respect to the arc length $s$.
In this paper, instead of directly following the Frenet frame definition, we
directly rotate the tangent vector $\alpha(s)$ counterclockwise for
$90^{\circ}$ in the plane of the axial curve so as to get the unit normal
vector at this point, i.e.,

$\bm{\beta}(s)=(-Y^{\prime}(s),X^{\prime}(s),0),$ (3)

and the binormal vector is accordingly defined as

$\bm{\gamma}(s)=(0,0,1)$ (4)

following right-hand rule. Now on the ground of this unit orthogonal frame, we
are able to develop a new orthogonal coordinate system in terms of $(s,y,z)$
and thereby the vector function of a specific point $(s,y,z)$ is

$\bm{r}(s,y,z)=\bm{r_{0}}(s)+y\bm{\beta}(s)+z\bm{\gamma}(s).$ (5)

###  II.2 Solution to the PB equation in terms of new curvilinear coordinates

We consider a 2D crooked nanochannel with charged wall containing binary
symmetric electrolyte solution(charge number $z_{\pm}=\pm 1$), as shown in
Fig. LABEL:fig1. The nanochannel is connected to two reservoirs with a fixed
distance apart as sources of ions (Fig. LABEL:fig1(b)). In order to study the
nanochannel curvature dependence, we let its contour length be variable.
According to the Gouy-Chapman’s model, counterions of the surface charges are
strongly adsorbed to the channel wall, forming a bounded layer of ions, and
the concentration of other free counterions tends to decline by the distance
from the interface so as to maintain electroneutrality. This screening region
with two zones, i.e., the zone of bounded layer and the zone of diffusion
layer, is EDL. In electrostatic situations, it is justifiable to consider the
ions in EDL to be in quasithermal equilibrium state when no biased potential
or extra pressure is applied to the nanochannel, hence the concentration of
both kinds of ions is supposed to obey the Boltzmann distribution as a
function of the electric potential, i.e.,

$c_{\pm}=c_{\rm b}\exp\bigg{(}-\frac{z_{\pm}e\psi}{k_{\rm B}T}\bigg{)},$ (6)

where $c_{b}$, $e$, $k_{\rm B}$, and $T$ are the molar bulk concentration far
away from the charged channel wall, electric charge, the Boltzmann constant
and the temperature (which is set to be 298.15 K in the present calculation)
respectively. In addition, on the ground of mean-field assumption, ions in
dilute electrolyte solution can be regarded as point charges with their
distribution determined by the surrounding electric field, making the
electrostatic equation a Poisson-type, which is given by

$\nabla^{2}\psi=-\frac{\rho}{\varepsilon_{0}\varepsilon_{\rm r}},$ (7)

where $\rho=F(c_{+}-c_{-})$, $\varepsilon_{0}$ and $\varepsilon_{\rm r}$ are
the volume charge density, the permittivity of vacuum and the relative
permittivity of the electrolyte respectively. Here $F$ is the Faraday
constant. Substituting Eq. (6) into Eq. (7) leads to the PB equation. However,
except for some 1D cases, it is extremely hard to solve the PB equation
analytically and practically a so-called Debye-Hückel (DH) approximation with
the premise of low surface electric potential ($\psi_{\rm s}<25$ mV) is
frequently used to linearize the PB equation. In DH approximation, expanding
the exponential term of electric potential and only keeping the first order
term leads to

$\nabla^{2}\psi=k^{2}\psi,$ (8)

where $k$, for binary electrolyte, is given by

$k=\sqrt{\frac{2F^{2}c_{\rm b}}{\varepsilon_{0}\varepsilon_{\rm r}RT}}.$ (9)

Here $k$ and $R$ are the DH parameter, also defined as the reciprocal of Debye
length $\lambda_{\rm D}$, which is the characteristic distance of EDL, and the
universal gas constant respectively.

In the curvilinear coordinate system that we develop in the previous section,
the scale factors, also called Lamé coefficients, are defined as

$h_{s}=1-yK(s),h_{y}=1,h_{z}=1,$ (10)

where $K(s)$ stands for the curvature of a planar curve defined in real number
field, i.e.,

$K(s)^{2}=X^{\prime\prime}(s)^{2}+Y^{\prime\prime}(s)^{2}.$ (11)

To be exact, $K(s)$ is positive when the normal vector $\bm{\beta}(s)$ at
point $s$ directs toward the center of curvature, and negative when
$\bm{\beta}(s)$ points to the opposite direction. With these coefficients, the
Laplacian for electric potential can be expressed in terms of these
curvilinear coordinates. The 2D crooked nanochannel we consider has a variable
contour length $L$ and a width $2a=8$ nm. Its central axis is chosen as an
asymmetric one defined by two parametric equations, in terms of $t$, as

$\left\\{\begin{aligned} X=&200t\\\ Y=&m\sin(\cfrac{5}{2}\pi
t^{2}),\end{aligned}\right.$ (12)

where $m$ is a parameter controlling its profile. In the 2D nanochannel placed
in $XY$ plane, the curvilinear coordinate system can be reduced to a 2D
fashion using space variables of arc length $s$ and transverse coordinate $y$.
Since typically contour length $L$ of a nanochannel is much larger than its
width $2a$, we let $L/a$ be over $50$ so that the curvature $K(s)$ is small if
no kinks exist. Therefore it is reasonable to consider only the first-order
terms of $yK$. For simplicity and without losing accuracy, the crooked shape
of the curved channel can be treated as a deviation from a straight channel,
which enables us to apply perturbation theory to solve the PB equation.

First of all we transform the space variables to a dimensionless form in
following rules as $\bar{x}=y/a,\bar{s}=s/L,\bar{K}=KL,\bar{k}=ka$ and define
a small quantity $\xi=a/L$ for perturbation calculation. Then the PB equation,
expanded to the first order of $\xi$, is given by

$\left(\frac{\partial^{2}}{\partial\bar{y}^{2}}-\bar{k}^{2}\right)\psi+\xi\left(-\bar{y}\bar{K}\frac{\partial^{2}}{\partial\bar{y}^{2}}-\bar{K}\frac{\partial}{\partial\bar{y}}+\bar{y}\bar{K}\bar{k}^{2}\right)\psi=0.$
(13)

Correspondingly, we expand the electric potential in terms of $\xi$, keeping
the prime and the first-order terms as

$\psi\approx\psi_{0}(\bar{y})+\xi\psi_{1}(\bar{y},\bar{s}).$ (14)

Equating the same order terms, we achieve the general solution to the PB
equation as

$\displaystyle\psi_{0}(\bar{y})=c_{1}e^{\lambda}+c_{2}e^{-\lambda}$ (15)
$\displaystyle\psi_{1}(\bar{y},\bar{s})=\frac{\bar{K}}{4\bar{k}}[(2c_{1}\lambda+c_{3})e^{\lambda}+(2c_{2}\lambda+c_{4})e^{-\lambda}],$
(16)

where $\lambda=\bar{k}\bar{y}$ and $c_{1},c_{2},c_{3},c_{4}$ are all constants
to be determined by boundary conditions.

##  III Numerical simulation

Strictly speaking, the general equations that govern electrokinetic effects in
a nanochannel as shown in Fig. LABEL:fig1(b) are Poisson-Nernst-Planck(PNP)
equations [16, 13, 14, 5]. As the nanochannel reaches the steady state with
little convection, the PNP equations combined with continuity equation become

$\displaystyle\nabla^{2}\psi=-\frac{F(c_{+}-c_{-})}{\epsilon}$ (17)
$\displaystyle\nabla\cdot\left(D_{\pm}\nabla
c_{\pm}+\frac{D_{\pm}Fz_{\pm}c_{\pm}}{RT}\nabla\psi\right)=0,$ (18)

where $D_{\pm}$ is the diffusion coefficients of cations and anions (for
simplicity, we assume they are constants[32, 33]). Note that the PNP equations
reduce to the PB equation in thermal equilibrium, a state without flux. The
coupled partial differential equations Eqs. (17) and (18) are numerically
solved by using finite element method program COMSOLTM. The two insulated
reservoirs connected

Figure 2: (a) Surface potential versus bulk concentration ($\psi$-$c_{\rm b}$)
for different average surface charge densities ($\sigma_{\rm a}=Q_{0}/S$).
Note that the vertical axis is in reverse order; (b) Surface potential versus
average surface charge density ($\psi$-$\sigma$). Line: analytical PB
solution; Dot: PNP solution via numerical calculation.

to the curved nanochannel, acting as the sources of electrolyte,
mathematically play the role of zero electric potential and constant
concentration boundary for Eqs. (17) and (18). During the process of
simulation, the two nanochannel walls are considered as charged surfaces in
different boundary conditions, and we fix the distance between the two
reservoirs and adjust the profile of the nanochannel via changing the
parameter $m$.

Figure 3: Surface charge density - arc length($\sigma$-s) distribution in
different bulk concentrations. Here we only show the surface charge density at
the bottom side of the nanochannel because that at the top is nothing but a
symmetric one about the black horizontal dashed line ($\sigma=-0.02~{}{\rm
e/nm^{2}}$). Inset: curvature-arc length (K-s) curve.

##  IV Results and discussion

###  IV.1 Conductive channel walls

We consider two kinds of boundary conditions for the PB equation,
corresponding to different mechanisms of the interaction between nanochannel
and electrolyte. Firstly, a conductive boundary, which can be found in
lamellar graphene nanochannels [25, 26], is considered for the 2D nanochannel.
An equipotential distribution of the two walls will maintain when they contact
with electrolyte solution. More specifically if we consider an asymmetrically
crooked nanochannel with each wall carrying a total amount of charge $Q_{0}$,
the boundary condition can be written as

$\int_{S}\vec{D}\cdot d\vec{A}=Q_{0},$ (19)

where $S$, $\vec{D}$ and $\vec{A}$ stand for the surface area, the electric
displacement, and the area vector respectively. Figure 2 shows how surface
potential of the channel depends on the average surface charge density and
concentration of ions in electrolyte solution, where both analytical (lines)
and numerical (dots) calculation results are illustrated. A closer look at
Fig. 2 exhibits a power law dependence between surface potential and bulk
concentration

$\psi|_{y=\pm a}\varpropto-c_{\rm b}^{-0.544},$ (20)

which is found to decrease linearly with the increase of the total amount of
surface charge (Fig. 2). Such a dependence on bulk concentration comes from
the screen effect of the EDL since higher ionic concentration means more ions
can be transferred to the surface so as to balance the excess charges. As
expected, most analytical results match well with the numerical simulation
results, except for the cases of low bulk concentrations and high surface
charge densities.

In order to investigate the electrostatic features along the nanochannel, we
also plot the charge density distribution against the arc length, as shown in
Fig. 3. As expected, there exists a tendency of point accumulation at the
places of big curvatures, which is found to enhance if the ionic concentration
decreases. When the ionic concentration increases, the screen effect of
counterions grows stronger and as a consequence the surface charge tends to
uniformly distribute along the channel, pushing the distribution curves closer
to the black horizontal dashed line, as shown in Fig. 3. Conclusively, the
profile of the nanochannel highly affects the EDL since the Debye length will
shorten with the increase of the bulk electrolyte concentration, weakening the
boundary effect.

Figure 4: Concentration of cations (above) and anions (below) as a function of
transverse coordinate($c_{\pm}-y$) at points of different curvatures. Black
dashed line indicates the bulk concentration of both kinds of ions.

###  IV.2 Uniformly Charged Channel Walls

In some biological and artificial systems, the channel wall is coated with
functional groups to control transportational properties of a nanochannel,
like permeability and permselectivity [34, 35, 36]. For example, nanochannels
made of silicate or boron nitride [37] after surface chemical modification
will carry immobile

Figure 5: Surface potential versus curvature($\psi$-$K$). The surface charge
density is set as -0.01 e/nm2 and the bulk concentration is 10 mM. Red and
blue lines correspond to top and bottom channel wall respectively. Solid line:
PB solution. Dotted line: PNP solution

surface charge in contact with electrolyte solution and theoretically the
boundary condition can be considered as fixed surface charge density, like
$\sigma_{0}$, as

$\left.\mp\varepsilon_{0}\varepsilon_{\rm r}\frac{\partial\psi}{\partial
y}\right|_{\pm a}=\sigma_{0}.$ (21)

Figure 4 shows the analytical result of transverse ionic concentration
distribution in the channel at places of different curvatures. As expected,
the concentration of cations, decaying from the two channel walls to its axis
as a typical feature of EDL, is higher than that of anions in negatively
charged nanochannel. Naturally at points of zero curvatures, ionic
concentration distribution is symmetric in $y$ direction, while at non-zero
curvature points, cations tend to accumulate near the surface that protrudes
outward. This preference for staying at concave areas of the channel
originates from profile-dependent surface electric potential. Even though the
charge density is uniform along the channel wall, the surface electric
potential actually varies with the curvature, inducing an asymmetric ionic
concentration distribution. Figure 5 presents the linear relationship between
surface electric potential and curvature, which provides a possibility of
controlling ion distribution in nanochannel via quantitatively adjusting its
profile. The reason why the discrepancy between the analytical and the
numerical results grow larger as the curvature increases is that high-order
terms of curvature in the PB equation is no longer negligible when the
curvature radius is reduced to the same scale as the width of the channel.

##  V Conclusion

We have developed a curvilinear coordinate system for a curved nanochannel,
which provides a straightforward and convenient mathematical representation of
the static and dynamic process inside. The analytical solutions to the 2D PB
equation quantify the profile effect on the electric double layer, which can
be seen by the dependence of both surface charge and ionic concentration
distributions on the curvature of the nanochannel. These results confirm the
possibility to dynamically control the ion distribution as well as the
electric potential through bending the nanochannel, a technique that may find
its promising application in flexible electronic devices. Further
investigation on the curvature effect in 3D nanofluidics can be made after the
mathematical challenge of solving the PB equation is overcome. Moreover, a
biased potential can be imposed to the two ends of a curved nanochannel so as
to induce a directional movement of ions, i.e., an electric current. It has
been reported that a rectification on the voltage can be seen in an
asymmetrically bent nanochannel [22, 23] and the mechanism behind the
phenomenon is yet to be clarified.

###### Acknowledgements.

This work was funded by the National Science Foundation of China under Grant
No. 11974292 and No. 11947401.

## References

  * [1] C. Duan and A. Majumdar, Nature nanotechnology 5, 848 (2010).
  * [2] X. Hou, W. Guo, and L. Jiang, Chemical Society Reviews 40, 2385 (2011).
  * [3] J. Geng, K. Kim, J. Zhang, A. Escalada, R. Tunuguntla, L. R. Comolli, F. I. Allen, A. V. Shnyrova, K. R. Cho, D. Munoz, et al., Nature 514, 612 (2014).
  * [4] R. J. Hunter, Foundations of colloid science (Oxford university press, 2001).
  * [5] R. B. Schoch, J. Han, and P. Renaud, Reviews of modern physics 80, 839 (2008).
  * [6] J. C. Berg, An introduction to interfaces & colloids: the bridge to nanoscience (World Scientific, 2010).
  * [7] M. Gouy, (1910).
  * [8] D. L. Chapman, The London, Edinburgh, and Dublin philosophical magazine and journal of science 25, 475 (1913).
  * [9] R. J. Hunter, Colloid science: zeta potential in colloid science: principles and applications (1981).
  * [10] C. L. Rice and R. Whitehead, The Journal of Physical Chemistry 69, 4017 (1965).
  * [11] H. Ohshima, T. W. Healy, and L. R. White, Journal of colloid and interface science 90, 17 (1982).
  * [12] H. Ohshima, Colloid and Polymer Science 296, 647 (2018).
  * [13] I. Vlassiouk, S. Smirnov, and Z. Siwy, Nano letters 8, 1978 (2008).
  * [14] I. Vlassiouk, S. Smirnov, and Z. Siwy, Acs Nano 2, 1589 (2008).
  * [15] S. Movahed and D. Li, Electrophoresis 32, 1259 (2011).
  * [16] D. Constantin and Z. S. Siwy, Physical Review E 76, 041202 (2007).
  * [17] Y. Green, Y. Edri, and G. Yossifon, Physical Review E 92, 033018 (2015).
  * [18] Y. Green, Physical Review E 98, 033114 (2018).
  * [19] A. Rustom, R. Saffrich, I. Markovic, P. Walther, and H.H. Gerdes, Science 303, 1007 (2004).
  * [20] N. M. Sherer, M. J. Lehmann, L. F. Jimenez-Soto, C. Horensavitz, M. Pypaert, and W. Mothes, Nature cell biology 9, 310 (2007).
  * [21] D. Huh, K. Mills, X. Zhu, M. A. Burns, M. Thouless, and S. Takayama, Nature materials 6, 424 (2007).
  * [22] M. Wang, H. Meng, D. Wang, Y. Yin, P. Stroeve, Y. Zhang, Z. Sheng, B. Chen, K. Zhan, and X. Hou, Advanced Materials 31, 1805130 (2019).
  * [23] M. Wang, Y. Hou, L. Yu, and X. Hou, Nano Letters (2020).
  * [24] M.-L. Liu, M. Huang, L.-Y. Tian, L.-H. Zhao, B. Ding, D.-B. Kong, Q.-H. Yang, and J.-J. Shao, ACS applied materials & interfaces 10, 44915 (2018).
  * [25] A. R. Koltonow and J. Huang, Science 351, 1395 (2016).
  * [26] C. Cheng, G. Jiang, C. J. Garvey, Y. Wang, G. P. Simon, J. Z. Liu, and D. Li, Science advances 2, e1501272 (2016).
  * [27] J. H. Yun, M.-S. Chun, and H. W. Jung, Physics of Fluids 22, 052004 (2010).
  * [28] M.-S. Chun, Physical Review E 83, 036312 (2011).
  * [29] D. J. Struik, Lectures on classical differential geometry (Courier Corporation, 1961).
  * [30] R. L. Ricca, Fluid Dynamics Research 36, 319 (2005).
  * [31] d. A. Garcia et al., Physics of Plasmas 13, 022309 (2006).
  * [32] M. S. Kilic, M. Z. Bazant, and A. Ajdari, Physical review E 75, 021502 (2007).
  * [33] M. S. Kilic, M. Z. Bazant, and A. Ajdari, Physical review E 75, 021503 (2007).
  * [34] D. Stein, M. Kruithof, and C. Dekker, Physical Review Letters 93, 035901 (2004).
  * [35] R. H. Tunuguntla, R. Y. Henley, Y.-C. Yao, T. A. Pham, M. Wanunu, and A. Noy, Science 357, 792 (2017).
  * [36] S. J. Kim, Y.-C. Wang, J. H. Lee, H. Jang, and J. Han, Physical review letters 99, 044501 (2007).
  * [37] A. Siria, P. Poncharal, A.-L. Biance, R. Fulcrand, X. Blase, S. T. Purcell, and L. Bocquet, Nature 494, 455 (2013).