2022

[1]E.I. Jafarov

[1]Institute of Physics, Azerbaijan National Academy of Sciences, Javid av.
131, Baku, AZ1143, Azerbaijan

# On the exactly-solvable semi-infinite quantum well of the non-rectangular
step-harmonic profile

<EMAIL_ADDRESS>S.M. Nagiyev<EMAIL_ADDRESS>*

###### Abstract

An exactly-solvable model of the non-relativistic harmonic oscillator with a
position-dependent effective mass is constructed. The model behaves itself as
a semi-infinite quantum well of the non-rectangular profile. Such a form of
the profile looks like a step-harmonic potential as a consequence of the
certain analytical dependence of the effective mass from the position and
semiconfinement parameter $a$. Both states of the discrete and continuous
spectrum are studied. In the case of the discrete spectrum, wavefunctions of
the oscillator model are expressed through the Bessel polynomials. The
discrete energy spectrum is non-equidistant and finite as a consequence of its
dependence on parameter $a$, too. In the case of the continuous spectrum,
wavefunctions of the oscillator model are expressed through the ${}_{1}F_{1}$
hypergeometric functions. At the limit, when the parameter $a$ goes to
infinity, both wavefunctions, and the energy spectrum of the model under
construction correctly reduce to corresponding results of the usual non-
relativistic harmonic oscillator with a constant effective mass. Namely,
wavefunctions of the discrete spectrum recover wavefunctions in terms of the
Hermite polynomials, and wavefunctions of the continuous spectrum simply
vanish. We also present a new limit relation that reduces Bessel polynomials
directly to Hermite polynomials and prove its correctness using the
mathematical induction technique.

###### keywords:

Semi-infinite potential, Exact solution, Harmonic oscillator, Position-
dependent effective mass, Bessel polynomials

##  1 Introduction

Low-dimensional quantum systems, especially quantum well systems are of great
interest in both theoretical and experimental physics as well as related areas
due to the recent impact of such quantum systems on the development of modern
nanotechnologies and advanced devices. Many physical phenomena related to
modern nanotechnologies can be successfully explained through the exact or
numerical solution of one-dimensional Schrödinger equation with effective mass
$m_{0}$ and time-independent potential $V\left(x\right)$. Quantum well systems
are among them, too. Most of the quantum well potentials leading to the exact
solution of the corresponding Schrödinger equation are symmetric in respect to
position inversion due to that, such symmetric behavior of the potential
substantially simplifies the exact solution of the Schrödinger equation for
it. However, the recent development of the molecular-beam epitaxy and related
pioneering methods on the experimental fabrication of the artificial quantum
well structures with profiles differing from the symmetric ones miller1984a ;
miller1984b ; miller1985 ; gossard1986  also requires to describe accurately a
number of experimental results in the nanotechnologies via the exactly-
solvable potentials behaving themselves as confined quantum well systems with
semi-infinite property and non-rectangular profile.

In the present paper, we report on the new quantum oscillator model that
behaves itself as a semi-infinite quantum well of the non-rectangular profile.
If to be exact, then this quantum system is confined by an infinitely high
wall at an arbitrary negative finite value of a position and with a finite
wall at an arbitrary positive finite value of the position. Also, the system
behaves itself as a non-relativistic quantum harmonic oscillator within these
two values of the position, therefore, it has a completely non-rectangular
profile. In other words, its behavior almost overlaps with the confined
quantum system under the so-called step-harmonic potential rizzi2010 ;
amthong2014 , but differs from it that the model under the current study does
not vary abruptly. Such a behavior of the quantum system is achieved thanks to
the introduction of the effective mass varying with position morris2015 ;
morris2017 .

We structured our paper as follows: the next section presents a brief review
of the exact solution of the non-relativistic quantum harmonic oscillator
within the canonical approach. Then, Sect. 3 consists of the main result of
the current paper – it is mainly devoted to the construction of the exactly
solvable model of the non-relativistic harmonic oscillator with a position-
dependent effective mass, which behaves itself as a semi-infinite quantum well
of the non-rectangular profile. In the final section, we discuss the limit
cases, when the parameter $a$ goes to infinity and as a consequence, the
wavefunctions of the model under construction expressed by the Bessel
polynomials completely recover the wavefunctions of the non-relativistic
quantum harmonic oscillator within the canonical approach as well as energy
spectrum of non-equidistant and finite form becomes equidistant and infinite.

##  2 The non-relativistic quantum harmonic oscillator within the canonical
approach

In this informative section, we are going to provide briefly mathematical
expressions belonging to the non-relativistic quantum harmonic oscillator
within the canonical approach. These expressions also can be found easily in
most quantum mechanics textbooks. But, we decided to include them in our
paper, because, they will be useful during the computations being done in the
next section.

Wavefunctions of the stationary states and the energy spectrum of the non-
relativistic quantum harmonic oscillator within the canonical approach can be
obtained exactly from the one-dimensional time-independent Schrödinger
equation

$\hat{H}\psi\left(x\right)=E\psi\left(x\right),$ (2.1)

where the full Hamiltonian is a sum of the kinetic and potential energy
operators

$\hat{H}=\hat{H}_{0}+V\left(x\right),$ (2.2)

with the kinetic energy operator of the following general form:

$\hat{H}_{0}=\frac{{\hat{p}_{x}^{2}}}{{2m_{0}}}.$ (2.3)

Here, the momentum operator $\hat{p}_{x}$ being defined within the canonical
approach is known as

$\hat{p}_{x}=-i\hbar\frac{\mathrm{d}}{{\mathrm{d}x}}.$ (2.4)

We are going to deal with the model based on the harmonic oscillator,
therefore, its potential $V\left(x\right)$ with states bounded at infinity is
defined as

$V\left(x\right)=\frac{{m_{0}\omega^{2}x^{2}}}{2},\quad-\infty<x<+\infty.$
(2.5)

Here, $\omega$ is its constant angular frequency. Substitution of (2.3)–(2.5)
at (2.2) leads to the following second order differential equation generated
from Eq.(2.1):

$\frac{{\hbar^{2}}}{{2m_{0}}}\frac{{\mathrm{d}^{2}\psi}}{{\mathrm{d}x^{2}}}+\left({E-\frac{{m_{0}\omega^{2}x^{2}}}{2}}\right)\psi=0.$
(2.6)

Exact solution of this equation allows to obtain equidistant energy spectrum
$E\equiv E_{n}$ and wavefunctions of the stationary states
$\psi\left(x\right)\equiv\psi_{n}\left(x\right)$ of the following analytical
form:

$E_{n}=\hbar\omega\left({n+\frac{1}{2}}\right),\quad n=0,1,2,\ldots,$ (2.7)

$\psi_{n}\left(x\right)=\frac{1}{{\sqrt{2^{n}n!}}}\left({\frac{{m_{0}\omega}}{{\pi\hbar}}}\right)^{{\textstyle{1\over
4}}}\mathrm{e}^{-\frac{{m_{0}\omega
x^{2}}}{{2\hbar}}}H_{n}\left({\sqrt{\frac{{m_{0}\omega}}{\hbar}}x}\right).$
(2.8)

Here, $H_{n}\left({x}\right)$ are Hermite polynomials defined in terms of the
${}_{2}F_{0}$ hypergeometric functions as follows koekoek2010 :

$H_{n}(x)=(2x)^{n}\,_{2}F_{0}\left({\begin{array}[]{*{20}c}{\begin{array}[]{*{20}c}{-n/2,-(n-1)/2}\\\
-\\\ \end{array}\mathchar 24635\relax}&{-\frac{1}{{x^{2}}}}\\\
\end{array}}\right).$ (2.9)

These polynomials also satisfy the following recurrence relations, which we
will use in the final section of the current paper:

$H_{n+1}\left(x\right)=2xH_{n}\left(x\right)-2nH_{n-1}\left(x\right).$ (2.10)

$\psi_{n}\left(x\right)$ wavefunctions (2.8) are orthonormalized and they
satisfy the following orthogonality relation:

$\int\limits_{-\infty}^{\infty}{{\psi_{m}^{*}(x)}{\psi}_{n}(x)dx}=\delta_{mn}.$
(2.11)

Hamiltonian (2.2) being sum of the non-relativistic canonical kinetic (2.3)
and harmonic potential energy operators (2.5) can be easily factorized in
terms of the harmonic oscillator creation and annihilation operators as
follows dirac1927 ; infeld1951 ; messiah1966 :

$\hat{H}=\hbar\omega\left({\hat{a}^{+}\hat{a}^{-}+\frac{1}{2}}\right).$ (2.12)

Here, harmonic oscillator creation and annihilation operators have the
following first-order differential analytical expression:

$\displaystyle\hat{a}^{+}=\frac{1}{{\sqrt{2}\lambda_{0}}}\left({\lambda_{0}^{2}x-\frac{\mathrm{d}}{{\mathrm{d}x}}}\right),$
(2.13)
$\displaystyle\hat{a}^{-}=\frac{1}{{\sqrt{2}\lambda_{0}}}\left({\lambda_{0}^{2}x+\frac{\mathrm{d}}{{\mathrm{d}x}}}\right),$

where, $\lambda_{0}=\sqrt{\frac{{m_{0}\omega}}{\hbar}}$.

We note here two main properties of these operators, one of which is that they
do not commute, therefore their commutation is one, i.e.

$\left[{\hat{a}^{-},\hat{a}^{+}}\right]=1,$

and the second one is the following equation, where the action of the
annihilation operator $\hat{a}^{-}$ to the ground state of wavefunction (2.8)
is zero:

$\hat{a}^{-}\psi_{0}\left(x\right)=0.$

##  3 The model of the semi-infinite quantum well of the non-rectangular
profile

In this section, we are going to construct the exactly solvable model of the
non-relativistic harmonic oscillator with a position-dependent effective mass,
which behaves itself as a semi-infinite quantum well of the non-rectangular
profile. In previous section, it was highlighted that the non-relativistic
harmonic oscillator with a potential defined via (2.5) has the wavefunctions
(2.8), which vanish at $\pm\infty$. Then, we decided that if we want to make
from the non-relativistic harmonic oscillator a semi-infinite quantum well of
the non-rectangular profile in the form of the step-harmonic potential, here
the effective mass certainly depending on the position can be a powerful tool.
Before showing a certain implementation of this idea, we want to share some
known information about the conception of the position-dependent effective
mass in quantum mechanics. This conception is presented in bendaniel1966 ,
which supposed that the change of the band structure of the independent-
particle model of the superconductor harrison1961  on the tunneling experiment
giaever1960a ; giaever1960b  should be simulated by a spatially varying
effective mass $M\left(x\right)$. In other words, it is necessary to implement
a replacement $m_{0}\to M\left(x\right)$ in full Hamiltonian (2.2) and then to
introduce a new analogue of the kinetic energy operator (2.3). The following
simpler realization of the Hermitian kinetic energy operator, now called as
BenDaniel-Duke kinetic energy operator has been introduced:

$\hat{H}_{0}\equiv\hat{H}_{0}^{BD}=-\frac{{\hbar^{2}}}{2}\frac{\mathrm{d}}{{\mathrm{d}x}}\frac{1}{{M(x)}}\frac{\mathrm{d}}{{\mathrm{d}x}}.$
(3.1)

Figure 1:  Behavior of the semi-infinite quantum well potential (3.3) of the
non-rectangular step-harmonic profile (solid line) and the effective mass
$M\left(x\right)$ varying by position (3.4) (dashed line) for the confinement
parameter: a) $a=1$; b) $a=2$; c) $a=3$; d) $a=4$ ($m_{0}=\omega=\hbar=1$,
hence, $\lambda_{0}=1$).

Despite that there is a number of other definitions for the position-dependent
effective mass kinetic energy operator gora1969 ; zhu1983 ; vonroos1983 ;
morrow1984 ; li1993 ; li2009 ; lima2012 ; nobre2015 ; mustafa2019 ;
elnabulsi2020 , we are going to use in our computations BenDaniel–Duke
definition (3.1) due to that from mathematics viewpoint, the Schrödinger
equation with the other kinetic energy operator definitions does not differ
significantly from BenDaniel–Duke one and the method used here can be later
easily applied for them, too. Also, it is shown in kolesnikov1999  that if to
modify BenDaniel–Duke kinetic energy operator (3.1), then additional singular
terms will appear in its modified version and they can lead to a discontinuity
of the wavefunction.

By substitution of (3.1) at (2.2) one has the following Schrödinger equation:

$\frac{{\hbar^{2}}}{2}\frac{\mathrm{d}}{{\mathrm{d}x}}\frac{1}{M}\frac{\mathrm{d}}{{\mathrm{d}x}}\psi+\left({E-V}\right)\psi=0.$
(3.2)

Here, by definition $M\equiv M\left(x\right)$ and $V\equiv V\left(x\right)$.
As a next step, one needs to introduce a potential $V\left(x\right)$ that will
behave itself as a semi-infinite quantum well of the non-rectangular profile.
In other words, the quantum system under study will be confined at an
arbitrary negative finite value of a position by an infinitely high wall
($x=-a$), but at an arbitrary positive finite value of position with a finite
wall and is non-rectangular profile between these two walls will be achieved
thanks to its behavior as a non-relativistic quantum harmonic oscillator
within these two values of the position. Therefore, it looks like a step-
harmonic potential, but with no abrupt change. Analytical expression of such a
potential $V\left(x\right)$ can be defined as follows:

$V\left(x\right)=\left\\{\begin{array}[]{ll}\frac{{m_{0}\omega^{2}a^{2}}}{2},&\hbox{
for }\quad x\to+\infty,\\\ \frac{{M\omega^{2}x^{2}}}{2},&\hbox{ for }\quad
x>-a,\\\ +\infty,&\hbox{ for }\quad x\leq-a,\\\
\end{array}\right.\qquad(a>0).$ (3.3)

It is clear that the depth of the semi-infinite well is as follows:

$V_{\infty}\equiv\mathop{\lim}\limits_{x\to+\infty}V\left(x\right)=\frac{{m_{0}\omega^{2}a^{2}}}{2}.$

Then, the basic principles of quantum mechanics require that the energy
spectrum is discrete if $E<V_{\infty}$ and it is continuous if $E>V_{\infty}$.
From definition of potential (3.3), it is also clear that the height of a
finite wall depends on the parameter $a$ that defines the negative value of
the position for which the potential becomes infinite. Such a general behavior
can be completely achieved via the following analytical definition of $M\equiv
M\left(x\right)$:

$M\left(x\right)=\frac{{a^{2}m_{0}}}{{\left({a+x}\right)^{2}}}.$ (3.4)

In fig.1, we present a comparative behavior of the semi-infinite quantum well
potential (3.3) of the non-rectangular step-harmonic profile and the effective
mass $M\left(x\right)$ varying by position (3.4) for different values of the
confinement parameter $a$, namely, for $a=1$; $2$; $3$ and $4$. One observes
that the behavior of both functions extracted from Eqs. (3.3) and (3.4)
completely proves our main goal. The potential (3.3) behaves itself as a semi-
infinite quantum well potential of the non-rectangular step-harmonic profile
and the position-dependent effective mass defined analytically via Eq.(3.4)
directly defines its behavior as confined by an infinitely high wall at an
arbitrary negative finite value of a position and with a finite wall at an
arbitrary positive finite value of the position as well as a non-relativistic
quantum harmonic oscillator within these two values of the position. Greater
values of the parameter $a$ exhibit evidence of the disappearance of the
confinement effect from both negative and positive sides of the quantum-well
system.

Substitution of (3.3) and (3.4) at (3.2) leads to the following second-order
differential equation:

$\frac{{\mathrm{d}^{2}\psi}}{{\mathrm{d}x^{2}}}+\frac{2}{{a+x}}\frac{{\mathrm{d}\psi}}{{\mathrm{d}x}}-\left({\frac{{\frac{{m_{0}^{2}\omega^{2}a^{4}}}{{\hbar^{2}}}x^{2}}}{{\left({a+x}\right)^{4}}}-\frac{{\frac{{2m_{0}a^{2}E}}{{\hbar^{2}}}}}{{\left({a+x}\right)^{2}}}}\right)\psi=0.$
(3.5)

One applies the transition to the dimensionless variable $\xi=\frac{x}{a}$
that slightly changes Eq.(3.5) to the following one:

$\frac{{\mathrm{d}^{2}\psi}}{{\mathrm{d}\xi^{2}}}+\frac{2}{{1+\xi}}\frac{{\mathrm{d}\psi}}{{\mathrm{d}\xi}}-\left({\frac{{\frac{{m_{0}^{2}\omega^{2}a^{4}}}{{\hbar^{2}}}\xi^{2}}}{{\left({1+\xi}\right)^{4}}}-\frac{{\frac{{2m_{0}a^{2}E}}{{\hbar^{2}}}}}{{\left({1+\xi}\right)^{2}}}}\right)\psi=0.$
(3.6)

Introduction of the following constants

$c_{0}=\frac{{2m_{0}a^{2}E}}{{\hbar^{2}}},\quad
c_{2}=c_{0}+\lambda_{0}^{4}a^{4},$ (3.7)

and some easy mathematical tricks slightly simplify (3.6) as follows:

$\psi^{\prime\prime}+\frac{2}{{1+\xi}}\psi^{\prime}+\frac{{c_{0}\left({1+\xi}\right)^{2}+\left({c_{0}-c_{2}}\right)\xi^{2}}}{{\left({1+\xi}\right)^{4}}}\psi=0.$
(3.8)

Here, $\psi^{\prime}\equiv\frac{{\mathrm{d}\psi}}{{\mathrm{d}\xi}}$ and
$\psi^{\prime\prime}\equiv\frac{{\mathrm{d}^{2}\psi}}{{\mathrm{d}\xi^{2}}}$.

We look for the solution of the wave function as follows:

$\psi=\varphi y,$ (3.9)

where, $\varphi\equiv\varphi\left(\xi\right)$ has the following general
analytical form:

$\varphi=\left({1+\xi}\right)^{A}\mathrm{e}^{\frac{B}{{1+\xi}}}.$ (3.10)

Here, one needs to find exactly both parameters $A$ and $B$ through the
substitution of (3.10) together with (3.9) at (3.8). Long straightforward
computations lead to the following statement that it is possible only if

$A=B=\varepsilon\lambda_{0}^{2}a^{2},\quad\varepsilon=\pm 1,$ (3.11)

thereof, one obtains for $\varphi$

$\varphi=\left({1+\xi}\right)^{\varepsilon\lambda_{0}^{2}a^{2}}\mathrm{e}^{\frac{{\varepsilon\lambda_{0}^{2}a^{2}}}{{1+\xi}}}.$
(3.12)

Now, one needs to define $\varepsilon$. We take into account that for our
quantum model the following limit relations have to hold:

$\mathop{\lim}\limits_{\xi\to+\infty}\varphi=\mathop{\lim}\limits_{\xi\mathop{\to}\limits-1^{+}}\varphi=0.$

This means that

$\varepsilon=-1.$

Then, the final analytical expression of $\varphi$ is as follows:

$\varphi=\left({1+\xi}\right)^{-\lambda_{0}^{2}a^{2}}\mathrm{e}^{\frac{{-\lambda_{0}^{2}a^{2}}}{{1+\xi}}}.$
(3.13)

Thanks to exact expression of (3.13), the second-order differential equation
for $\psi$ (3.8) transfers to the following second-order differential equation
for $y$:

$\left({\xi+1}\right)^{2}y^{\prime\prime}+2\left[{\left({1-\lambda_{0}^{2}a^{2}}\right)\xi+1}\right]y^{\prime}=\left({\lambda_{0}^{2}a^{2}-c_{0}}\right)y.$
(3.14)

Further, one needs to solve Eq.(3.14) to obtain exact eigenvalues of the
energy spectrum. First, it is necessary to take into account that the
polynomial $\sigma\left(\xi\right)=\left({\xi+1}\right)^{2}$ has multiple
roots. Therefore, following the method from nikiforov1988 , it is better to
apply the transition to the new variable $\tau$ as follows:

$\tau=\frac{1}{{\xi+1}}.$ (3.15)

Substitution of (3.15) at Eq.(3.14) leads to

$\frac{{\mathrm{d}^{2}y}}{{\mathrm{d}\tau^{2}}}+2\lambda_{0}^{2}a^{2}\frac{{1-\tau}}{\tau}\frac{{\mathrm{d}y}}{{\mathrm{d}\tau}}+\frac{{c_{0}-\lambda_{0}^{2}a^{2}}}{{\tau^{2}}}y=0.$
(3.16)

We look for its solution as follows:

$y=\tau^{\gamma}v\left(\tau\right),$ (3.17)

where, $\gamma$ is an arbitrary parameter, which will be defined below. As a
result of the simple computations, we obtain the following equation for
$v\left(\tau\right)$:

$\frac{{\mathrm{d}^{2}v}}{{\mathrm{d}\tau^{2}}}+\frac{{\mu-2\lambda_{0}^{2}a^{2}\tau}}{\tau}\frac{{\mathrm{d}v}}{{\mathrm{d}\tau}}+\frac{{\rho-\delta\tau}}{{\tau^{2}}}v=0,$
(3.18)

where $\mu=2\left({\gamma+\lambda_{0}^{2}a^{2}}\right)$,
$\rho=c_{0}-\lambda_{0}^{2}a^{2}+\gamma\left({\gamma-1}\right)+2\gamma\lambda_{0}^{2}a^{2}$
and $\delta=2\gamma\lambda_{0}^{2}a^{2}$.

Now we choose the constant parameter $\gamma$ in such a way that it satisfies
the condition $\rho=0$. Hence we get

$\gamma=\frac{{1-2\lambda_{0}^{2}a^{2}+e\sqrt{1+4\lambda_{0}^{4}a^{4}-4c_{0}}}}{2},\quad
e=\pm 1.$ (3.19)

As a consequence of this condition, one obtains from Eq.(3.18) the following
equation for $v\left(\tau\right)$:

$\tau\frac{{\mathrm{d}^{2}v}}{{\mathrm{d}\tau^{2}}}+\left({\mu-2\lambda_{0}^{2}a^{2}\tau}\right)\frac{{\mathrm{d}v}}{{\mathrm{d}\tau}}-\delta
v=0.$ (3.20)

Now, if to change the variable from $\tau$ to $z=2\lambda_{0}^{2}a^{2}\tau$,
then the equation above will overlap with the following hypergeometric type
equation landau1991 :

$z\frac{{\mathrm{d}^{2}v}}{{\mathrm{d}z^{2}}}+\left({\mu-z}\right)\frac{{\mathrm{d}v}}{{\mathrm{d}z}}-\gamma
v=0.$ (3.21)

Analytical solution to this equation in terms of the ${}_{1}F_{1}$
hypergeometric functions is well known:

$v=C_{1}\cdot\,_{1}F_{1}\left({\begin{array}[]{*{20}c}\gamma\\\ \mu\\\
\end{array}\mathchar
24635\relax\;2\lambda_{0}^{2}a^{2}\tau}\right)+C_{2}\cdot\tau^{1-\mu}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}{\gamma-\mu+1}\\\
{2-\mu}\\\ \end{array}\mathchar
24635\relax\;2\lambda_{0}^{2}a^{2}\tau}\right).$ (3.22)

###  3.1 The states of the discrete spectrum ($E<V_{\infty}$)

Condition $E<V_{\infty}$ requires for $\gamma$ parameter from definition
(3.19) that $e=+1$ is only possible case. Therefore, one simplifies

$\gamma=\frac{{1-2\lambda_{0}^{2}a^{2}+\sqrt{1+4\lambda_{0}^{4}a^{4}-4c_{0}}}}{2}.$
(3.23)

Moreover, the condition for the positivity of the root expression should be
satisfied, too:

$1+4\lambda_{0}^{2}a^{2}-4c_{0}=\left({1-2\lambda_{0}^{2}a^{2}}\right)^{2}+4\lambda_{0}^{2}a^{2}-4c_{0}>0,$
(3.24)

where, $\lambda_{0}^{2}a^{2}-c_{0}<0$. Due to that
$\mu=1+\sqrt{1+4\lambda_{0}^{4}a^{4}-4c_{0}}>1$, hence the finiteness
requirement at value $\tau=0$ for function $v\left(\tau\right)$ defined
through (3.22) leads to $C_{2}=0$. Then, we obtain an analytical expression of
function $y$ (3.17) as follows:

$y\left(\tau\right)=\tau^{\gamma}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}\gamma\\\
\mu\\\ \end{array}\mathchar 24635\relax\;2\lambda_{0}^{2}a^{2}\tau}\right).$

Analytical expression of the wavefunction (3.9) also can be written down as
follows:

$\psi\left(\xi\right)=C\left({1+\xi}\right)^{-\gamma-\lambda_{0}^{2}a^{2}}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{2}}}{{1+\xi}}}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}\gamma\\\
\mu\\\ \end{array}\mathchar
24635\relax\;2\frac{{\lambda_{0}^{2}a^{2}}}{{1+\xi}}}\right),\qquad C\equiv
C_{1}.$ (3.25)

It is finite at $\xi\to\infty$. Then, another requirement
$\psi\left(-1\right)=0$ leads to quantization of the energy through
$\gamma=-n$, i.e.

$\sqrt{1+4\lambda_{0}^{4}a^{4}-4c_{0}}=2\lambda_{0}^{2}a^{2}-2n-1,\quad
n=0,1,2,3,\ldots.$

Positivity condition (3.24) restricts $n$ from above, i.e. due to condition
$2\lambda_{0}^{2}a^{2}-2n-1>0$, there is a maximum value of $n$ as
$n_{max}=N=0,1,2,3,\ldots$. Here, $N$ itself satisfies the condition

$N<\lambda_{0}^{2}a^{2}-\frac{1}{2}.$ (3.26)

Now, taking into account above listed restrictions one easily obtains for the
energy spectrum that

$E\equiv
E_{n}^{QW}=\hbar\omega\left({n+\frac{1}{2}}\right)-\frac{{\hbar^{2}}}{{2m_{0}a^{2}}}n\left({n+1}\right),\quad
n=0,1,2,\ldots,N.$ (3.27)

Hence, an analytical expression of the wavefunctions of the stationary states
(3.25) also slightly changes as follows:

$\psi_{n}\left(\xi\right)=C_{n}\left({1+\xi}\right)^{n-\lambda_{0}^{2}a^{2}}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{2}}}{{1+\xi}}}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}{-n}\\\
{2\lambda_{0}^{2}a^{2}-2n}\\\ \end{array}\mathchar
24635\relax\;2\frac{{\lambda_{0}^{2}a^{2}}}{{1+\xi}}}\right),\quad
n=0,1,2,\ldots,N.$ (3.28)

Using the following well-known expression of the generalized Laguerre
polynomials through ${}_{1}F_{1}$ hypergeometric functions koekoek2010

$L_{n}^{\alpha}\left(x\right)=\frac{{\left({\alpha+1}\right)_{n}}}{{n!}}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}{-n}\\\
{\alpha+1}\\\ \end{array}\mathchar 24635\relax\;x}\right),$

analytical expression of the wavefunctions of the stationary states (3.28)
also can be written down through these polynomials as follows:

$\psi_{n}\left(\xi\right)=C_{n}\left({1+\xi}\right)^{n-\lambda_{0}^{2}a^{2}}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{2}}}{{1+\xi}}}L_{n}^{2\lambda_{0}^{2}a^{2}-2n-1}\left({2\frac{{\lambda_{0}^{2}a^{2}}}{{1+\xi}}}\right).$
(3.29)

One needs to note that there exist Bessel polynomials from the Askey scheme of
the orthogonal polynomials koekoek2010 , which are expressed through both
${}_{1}F_{1}$ and ${}_{2}F_{0}$ hypergeometric functions as follows:

$\displaystyle y_{n}\left(x\mathchar 24635\relax\;\alpha\right)$
$\displaystyle=$
$\displaystyle\left({n+\alpha+1}\right)_{n}\left({\frac{x}{2}}\right)^{n}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}{-n}\\\
{-2n-\alpha}\\\ \end{array}\mathchar 24635\relax\;\frac{2}{x}}\right)$ (3.32)
$\displaystyle=$
$\,{}_{2}F_{0}\left({\begin{array}[]{*{20}c}{-n,n+\alpha+1}\\\ -\\\
\end{array}\mathchar 24635\relax\;-\frac{x}{2}}\right),\quad
n=0,1,2,\ldots,N.$ (3.35)

Bessel polynomials $y_{n}\left({\frac{{2\xi+2d}}{\beta}\mathchar
24635\relax\;\alpha}\right)$ satisfy the second-order differential equation of
the following form koekoek2010

$\left({\xi+d}\right)^{2}y_{n}^{\prime\prime}\left(\xi\right)+\left[{\left({\alpha+2}\right)\left({\xi+d}\right)+\beta}\right]y_{n}^{\prime}\left(\xi\right)=n\left({n+\alpha+1}\right)y_{n}\left(\xi\right),$
(3.36)

where

$-d<\xi<+\infty,\quad\alpha<-\left({2N+1}\right),\quad\beta>0,$

and

$y_{n}\left(\xi\right)\equiv y_{n}\left({\frac{{2\xi+2d}}{\beta}\mathchar
24635\relax\;\alpha}\right).$

Therefore, wavefunctions (3.28) can also be written down in terms of the
Bessel polynomials, and the normalization factor $C_{n}$ for these
wavefunctions can be easily computed from the orthogonality relation of these
polynomials.

Taking into account that the following orthogonality relation holds for the
Bessel polynomials $y_{n}\left({\frac{{2\xi+2d}}{\beta}\mathchar
24635\relax\;\alpha}\right)$ koekoek2010 :

$\int\limits_{-d}^{+\infty}{\left({\xi+d}\right)^{\alpha}\mathrm{e}^{-\frac{\beta}{{\xi+d}}}y_{m}\left({\frac{{2\xi+2d}}{\beta}\mathchar
24635\relax\;\alpha}\right)y_{n}\left({\frac{{2\xi+2d}}{\beta}\mathchar
24635\relax\;\alpha}\right)\mathrm{d}\xi}=-\frac{{\beta^{\alpha+1}}{\Gamma\left({-n-\alpha}\right)n!}}{{2n+\alpha+1}}\delta_{mn},$

one finds the orthonormalized wavefunctions (3.28) or (3.29) in terms of the
Bessel polynomials, which are as follows:

$\psi_{n}^{QW}\left(x\right)=C_{n}^{QW}\left({\frac{x}{a}+1}\right)^{-\lambda_{0}^{2}a^{2}}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{3}}}{{x+a}}}y_{n}\left({\frac{{x+a}}{{\lambda_{0}^{2}a^{3}}}\mathchar
24635\relax\;-2\lambda_{0}^{2}a^{2}}\right),$ (3.37)

where the orthonormalized coefficients $C_{n}^{QW}$ are

$C_{n}^{QW}=\left({2\lambda_{0}^{2}a^{2}}\right)^{\lambda_{0}^{2}a^{2}}\sqrt{\frac{{2\lambda_{0}^{2}a^{2}-2n-1}}{{2\lambda_{0}^{2}a^{3}n!\Gamma\left({2\lambda_{0}^{2}a^{2}-n}\right)}}}.$
(3.38)

###  3.2 The states of the continuous spectrum ($E>V_{\infty}$)

The spectrum of the energy values, satisfying the condition $E>V_{\infty}$ is
continuous and spreads from $V_{\infty}=\frac{{m_{0}\omega^{2}a^{2}}}{2}$ to
$\infty$. In this case we have $c_{0}>\lambda_{0}^{4}a^{4}+\frac{1}{4}$.
Hence, the value of $\gamma$ defined via Eq.(3.19) and consequently, values of
$\mu$, $\rho$ and $\delta$ are complex. For example,

$\gamma=\frac{{1-2\lambda_{0}^{2}a^{2}+i\sqrt{4c_{0}-4\lambda_{0}^{4}a^{4}-1}}}{2}.$
(3.39)

Then, the eigenfunctions of the continuous spectrum will have the following
expression:

$\psi_{E}\left(x\right)={\mathop{\rm
const}\nolimits}\cdot\left({\frac{x}{a}+1}\right)^{-\gamma-\lambda_{0}^{2}a^{2}}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{3}}}{{x+a}}}\,_{1}F_{1}\left({\begin{array}[]{*{20}c}\gamma\\\
\mu\\\ \end{array}\mathchar
24635\relax\;\frac{{2\lambda_{0}^{2}a^{3}}}{{x+a}}}\right).$ (3.40)

By obtaining exact expressions for the discrete energy spectrum (3.27) and the
wavefunctions of the stationary states (3.37) and (3.40) we achieved our main
goal. We solved exactly the Schrödinger equation (3.2) corresponding to the
potential (3.3). We are going to discuss the main properties of the discrete
energy spectrum (3.27) and the wavefunctions of the stationary states (3.37)
corresponding to this discrete spectrum in the final section, where we will
also show that how both the energy spectrum (3.27) and the wavefunctions of
the stationary states (3.37) reduce to the energy spectrum (2.7) and the
wavefunctions of the stationary states (2.8) under the limit $a\to+\infty$.

##  4 Limit cases and conclusions

Recently, jafarov2021 ; jafarov2022  developed the exactly-solvable model of a
one-dimensional nonrelativistic canonical semiconfined quantum harmonic
oscillator with a mass that varies with position. Analytical definition of the
mass of that model differs from the analytical definition (3.4) with the power
of the denominator $(a+x)$ that equals to $1$. After, quesne2022  generalized
the same model to the case when the power of the denominator $(a+x)$ of the
position-dependent effective mass is greater than $0$, but, less than $2$. One
deduces from these results that if the power of the denominator $(a+x)$ with
less than $2$ still allows the quantum system to exhibit only the
semiconfinement effect, but, with the power of the denominator that equals to
$2$ behavior of the quantum system under the potential (3.3) changes
drastically and both infinitely high wall for the negative value and finite
wall for the positive value of the position appear.

Let’s discuss the results obtained within this paper in more detail. This will
allow us to understand the main differences of both the energy spectrum (3.27)
and the wavefunctions of the stationary states (3.37) from to the energy
spectrum (2.7) and the wavefunctions of the stationary states (2.8). First of
all, one needs to highlight the main difference – it arises during the
comparison of (3.14) with the equation for the Bessel polynomials (3.36) as
well as through condition (3.26), from where one obtains that

$N<\lambda_{0}^{2}a^{2}-\frac{1}{2}\Rightarrow
a>\frac{1}{{\lambda_{0}}}\sqrt{N+\frac{1}{2}}.$ (4.1)

One observes from the expression of the energy spectrum (3.27) that it is
discrete, but, unlike the energy spectrum (2.7), it is non-linear and finite,
where a number of its finite levels is defined by the semiconfinement
parameter $a$ via inequality (4.1). The origin of such behavior is due to that
if our potential is confined from the negative side by the infinitely high
wall, then the maximal height of the wall from the right side is
$\frac{{m_{0}\omega^{2}a^{2}}}{2}$ due to definition of the potential (3.3).
Due to that $\min\left(N\right)=0$, one observes that

$\min\left(a\right)>\frac{1}{{\sqrt{2}\lambda_{0}}}.$ (4.2)

This condition allows us to say that our definition $a>0$ introduced at the
initial stage during the definition of the potential (3.3) is not completely
accurate, therefore, should be replaced by a more accurate definition (4.2).

Figure 2: The semi-infinite quantum well potential (3.3) of the non-
rectangular step-harmonic profile, the corresponding energy levels (3.27) and
the probability densities $\lvert{\psi_{n}^{QW}(x)}\rvert^{2}$ of the
wavefunctions of the stationary states (3.37) for the ground and $N$ excited
states for the confinement parameter: a) $a=1$, so $N=0$; b) $a=2$, so $N=3$;
c) $a=3$, so $N=8$ and d) $a=4$, so $N=15$ ($m_{0}=\omega=\hbar=1$, hence,
$\lambda_{0}=1$).

In Fig.2, we depicted the semi-infinite quantum well potential (3.3) of the
non-rectangular step-harmonic profile, the corresponding energy levels (3.27)
and the probability densities $\lvert{\psi_{n}^{QW}(x)}\rvert^{2}$ of the
wavefunctions of the stationary states (3.37) for the ground and $N$ excited
states, where these states correspond for different values of the confinement
parameter $a$, namely, there are plots for $a=1$; $2$; $3$ and $4$. For
simplicity, all computations for these plots made for values
$m_{0}=\omega=\hbar=1$ ($\lambda_{0}=1$). One observes when $a=1$, that is a
special case of inequality (4.2) and corresponds to the plot 1(a), only the
ground state exists for the quantum system under construction. This is due to
that $\max\left(N\right)=0$ if $a=1$. Such a property of the quantum system
under consideration is very attractive due to that from a mathematical
viewpoint, the polynomial itself is one if $N=0$. Then, from a physics
viewpoint there is only the ground state energy level that equals to
$\frac{\hbar\omega}{2}$ and completely overlaps with the ground state energy
level of the non-relativistic quantum harmonic oscillator, but, the
wavefunction of the ground state differs from the wavefunction of the ground
state of the non-relativistic quantum harmonic oscillator that exhibits the
Gaussian distribution. Therefore, one observes from Fig.1a that the
probability of finding the quantum system close the infinitely high wall is
higher than at values $x=0$ and $x\to\frac{m_{0}\omega^{2}a^{2}}{2}$. As a
consequence of such a property, excited states appear and their number
increases coherently only with the increase of the value of the confinement
parameter $a$. Such a behavior can be easily observed from the plots (b)–(d)
of Fig. 2. We did not present here the plot when $a\to+\infty$ completely
recovers the so-called Hermite oscillator model, but, one can easily compute
the limit from the energy (3.27) to the energy spectrum (2.7). Also, it is
possible to show that there is a correct limit from the wavefunctions of the
stationary states (3.37) to the wavefunctions of the stationary states (2.8)
under the case $a\to+\infty$.

During computation of the above-mentioned limit relations, we are going to use
the following approximate ($\lvert x\rvert<<1$) and asymptotic ($\lvert
x\rvert\to\infty$) formulae:

$\displaystyle\sqrt{1+x}\approx 1+\frac{1}{2}x-\frac{1}{8}x^{2},$ (4.3)
$\displaystyle\frac{1}{{1+x}}\approx 1-x+x^{2},$ (4.4)
$\displaystyle\ln\left({1+x}\right)\approx x-\frac{1}{2}x^{2},$ (4.5)
$\displaystyle\Gamma\left(x\right)\cong\sqrt{\frac{{2\pi}}{x}}\mathrm{e}^{x\ln
x-x}.$ (4.6)

As we noted, one can easily show that the following limit relation is correct:

$\mathop{\lim}\limits_{a\to\infty}E_{n}^{QW}=E_{n}=\hbar\omega\left({n+\frac{1}{2}}\right).$

Here, one needs to highlight that ground state levels ($n=0$) of both energy
spectrum expressions (3.27) and (2.7) already overlap, i.e.:

$E_{0}^{QW}=E_{0}=\frac{1}{2}\hbar\omega.$

Let’s rewrite the wavefunctions of the stationary states (3.37) as follows:

$\psi_{n}^{QW}\left(x\right)=C_{n}^{QW}\varphi^{QW}\left(x\right)y_{n}^{QW}\left(x\right),$

where one has ($b=\lambda_{0}a$ is introduced to make simpler the expression
of the wavefunction)

$C_{n}^{QW}=\left({2b^{2}}\right)^{b^{2}}\sqrt{\frac{{\lambda_{0}\left({2b^{2}-2n-1}\right)}}{{2b^{3}n!\Gamma\left({2b^{2}-n}\right)}}},$
$\varphi^{QW}\left(x\right)=\left({\frac{{\lambda_{0}x}}{b}+1}\right)^{-b^{2}}e^{-\frac{{b^{3}}}{{\lambda_{0}x+b}}},$
$y_{n}^{QW}\left(x\right)=y_{n}\left({\frac{{\lambda_{0}x+b}}{{b^{3}}}\mathchar
24635\relax\;-2b^{2}}\right).$

Taking into account (4.3) one can rewrite that

$C_{n}^{QW}\approx\sqrt{\frac{{\lambda_{0}}}{{2bn!}}}\frac{{\left({2b^{2}}\right)^{b^{2}}}}{{\sqrt{\Gamma\left({2b^{2}-n}\right)}}}.$

Here, taking into account Stirling’s approximation (4.6) one has

$C_{n}^{QW}\approx\sqrt{\frac{{2^{n}\lambda_{0}}}{{n!\sqrt{\pi}}}}b^{n}\mathrm{e}^{b^{2}}.$
(4.7)

Now, let’s explore the asymptotical behavior of $\varphi^{QW}\left(x\right)$.
Using (4.4) and (4.5) and as a result of the straightforward computations one
obtains that

$\varphi^{QW}\left(x\right)=\mathrm{e}^{-b^{2}\ln\left({1+\frac{{\lambda_{0}x}}{b}}\right)}\mathrm{e}^{-b^{2}\left({1+\frac{{\lambda_{0}x}}{b}}\right)^{-1}}\cong\mathrm{e}^{b^{2}-\frac{1}{2}\lambda_{0}^{2}x^{2}}.$
(4.8)

To obtain the above expression, we did not take into account the terms of the
expansions, which do not give any contribution to the limit relation
$\mathop{\lim}\limits_{a\to\infty}\psi_{n}^{QW}\left(x\right)=\psi_{n}\left(x\right)$.

Combining together Eqs.(4.7) and (4.8) one obtains

$C_{n}^{QW}\varphi^{QW}\left(x\right)\cong\sqrt{\frac{{\lambda_{0}}}{{n!\sqrt{\pi}}}}\left({\sqrt{2}b}\right)^{n}\mathrm{e}^{-\frac{1}{2}\lambda_{0}^{2}x^{2}}.$
(4.9)

The correct limit that reduces the Bessel polynomials directly to the Hermite
polynomials is presented below as a proposition.

######  Proposition 1 (Limit relation from $y_{n}\left({x\mathchar
24635\relax\;\alpha}\right)$ to $H_{n}\left(x\right)$).

The Hermite polynomials $H_{n}\left(x\right)$ follow from the Bessel
polynomials $y_{n}\left({x\mathchar 24635\relax\;\alpha}\right)$ by setting
$x\to 2/\nu+2/\nu\sqrt{2/\nu}x$ and $\alpha\to-\nu$ and then letting
$\nu\to\infty$ in the following way:

$\mathop{\lim}\limits_{\nu\to\infty}\left({-1}\right)^{n}\left({2\nu}\right)^{n/2}y_{n}\left({\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x\mathchar
24635\relax\;-\nu}\right)=H_{n}\left(x\right).$ (4.10)

Proof: We prove the correctness of the above proposition using the
mathematical induction technique. We start from the following Rodrigues-type
formula for the Bessel polynomials koekoek2010 :

$y_{n}\left({x\mathchar
24635\relax\;\alpha}\right)=2^{-n}x^{-\alpha}\mathrm{e}^{\frac{2}{x}}\frac{{\mathrm{d}^{n}}}{{\mathrm{d}x^{n}}}\left({x^{2n+\alpha}\mathrm{e}^{-\frac{2}{x}}}\right).$

Then, one easily finds that

$y_{1}\left({x\mathchar
24635\relax\;\alpha}\right)=\frac{1}{2}\left({2+\alpha}\right)x+1,\quad
y_{2}\left({x\mathchar
24635\relax\;\alpha}\right)=\frac{1}{4}\left({3+\alpha}\right)\left({4+\alpha}\right)x^{2}+\left({3+\alpha}\right)x+1.$
(4.11)

Now, one can check that the limit relation (4.10) is correct for both
$y_{1}\left({x\mathchar 24635\relax\;\alpha}\right)$ and
$y_{2}\left({x\mathchar 24635\relax\;\alpha}\right)$, i.e.

$-\mathop{\lim}\limits_{\nu\to\infty}\sqrt{2\nu}\cdot
y_{1}\left({\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x\mathchar
24635\relax\;-\nu}\right)=2x=H_{1}\left(x\right),$
$\mathop{\lim}\limits_{\nu\to\infty}2\nu\cdot
y_{2}\left({\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x\mathchar
24635\relax\;-\nu}\right)=4x^{2}-2=H_{2}\left(x\right).$

Let’s prove the correctness of the limit relation (4.10) for case $n>2$. One
uses the following recurrence relation for the Bessel polynomials koekoek2010
:

$y_{n+1}\left({x\mathchar
24635\relax\;\alpha}\right)=A_{n}y_{n}\left({x\mathchar
24635\relax\;\alpha}\right)+B_{n}y_{n-1}\left({x\mathchar
24635\relax\;\alpha}\right),$ (4.12)

where,

$A_{n}=\frac{{\left({2n+\alpha+1}\right)\left[{2\alpha+\left({2n+\alpha}\right)\left({2n+\alpha+2}\right)x}\right]}}{{2\left({n+\alpha+1}\right)\left({2n+\alpha}\right)}},\;B_{n}=\frac{{n\left({2n+\alpha+2}\right)}}{{\left({n+\alpha+1}\right)\left({2n+\alpha}\right)}}.$
(4.13)

Now, let the limit relation (4.10) is correct for the polynomials
$y_{n}\left({z\mathchar 24635\relax\;-\nu}\right)$ and
$y_{n-1}\left({z\mathchar 24635\relax\;-\nu}\right)$ for arbitrary $n$ and
$z=\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x$. Then, this relation also
will be correct for case of the polynomial $y_{n+1}\left({z\mathchar
24635\relax\;-\nu}\right)$. To prove the correctness of this statement, one
needs to apply the substitutions
$z=\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x$ and $\alpha=-\nu$ at the
recurrence relation (4.12), next, multiply its both sides to the factor
$\left(-1\right)^{n+1}\left(2\nu\right)^{\frac{n+1}{2}}$ and then, take a
limit $\nu\to\infty$ from it. We will have

$\mathop{\lim}\limits_{\nu\to\infty}\left({-1}\right)^{n+1}\left({2\nu}\right)^{\frac{{n+1}}{2}}y_{n+1}\left({\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x\mathchar
24635\relax\;-\nu}\right)=\bar{A}_{n}H_{n}\left(x\right)+\bar{B}_{n}H_{n-1}\left(x\right).$
(4.14)

Here,

$\bar{A}_{n}=-\mathop{\lim}\limits_{\nu\to\infty}\sqrt{2\nu}A_{n},\;\bar{B}_{n}=\mathop{\lim}\limits_{\nu\to\infty}2\nu
B_{n}.$ (4.15)

By substituting (4.13) at (4.15) one easily computes that

$\bar{A}_{n}=2x,\;\bar{B}_{n}=-2n.$

Then, the limit relation (4.14) will be rewritten as follows:

$\mathop{\lim}\limits_{\nu\to\infty}\left({-1}\right)^{n+1}\left({2\nu}\right)^{\frac{{n+1}}{2}}y_{n+1}\left({\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x\mathchar
24635\relax\;-\nu}\right)=2xH_{n}\left(x\right)-2nH_{n-1}\left(x\right).$
(4.16)

Now, thanks to the recurrence relations for the Hermite polynomials (2.10),
one observes that

$\mathop{\lim}\limits_{\nu\to\infty}\left({-1}\right)^{n+1}\left({2\nu}\right)^{\frac{{n+1}}{2}}y_{n+1}\left({\frac{2}{\nu}+\frac{2}{\nu}\sqrt{\frac{2}{\nu}}x\mathchar
24635\relax\;-\nu}\right)=H_{n+1}\left(x\right).$

This proves the proposition on the correctness of the limit relation (4.10).

One needs to note that lesky1998  studies various properties of the Bessel
polynomials. The correct limit that reduces the Bessel polynomials directly to
the Hermite polynomials is among these studied properties. However, the limit
relation studied there was not suitable for application to the model under
present construction. Therefore, we decided to introduce a new limit relation
that reduces the Bessel polynomials directly to the Hermite polynomials
through (4.10) and it is exactly what we need to prove the correct reduction
of the wavefunction (3.37) to the wavefunction (2.8) if $a\to\infty$. Despite
that Bessel polynomials belong to the family of the polynomials introduced
within the finiteness property of the Jacobi-like polynomials routh1885 ;
romanovski1929  and there are correct limit relations from the pseudo-Jacobi
to Bessel polynomials koekoek2010  as well as from the pseudo-Jacobi to
Hermite polynomials jafarov2021b ; nagiyev2022 , but limit relation that
reduces the Bessel polynomials directly to the Hermite polynomials through
(4.10) completes all possible relations between these polynomials and it is
also possible to maintain the direct relation between it and the limit from
lesky1998 .

Now, combining together (4.9) and limit relation (4.10) between the Bessel and
Hermite polynomials, one observes that the wavefunction (3.37) also correctly
reduces to the wavefunction (2.8) if $a\to\infty$, i.e.

$\mathop{\lim}\limits_{a\to\infty}\psi_{n}^{QW}\left(x\right)=\psi_{n}\left(x\right).$

The wavefunctions of the continuous spectrum (3.40) simply vanish under the
limit $a\to\infty$. One can slightly change Eq.(3.39) in terms of $\mu$ as
follows:

$\mu=2\left({\gamma+\lambda_{0}^{2}a^{2}}\right)=1+i\sqrt{4c_{0}-4\lambda_{0}^{4}a^{4}-1}.$
(4.17)

Next, the condition $c_{0}>\lambda_{0}^{4}a^{4}+\frac{1}{4}$ leads to
statement that

$4c_{0}-4\lambda_{0}^{4}a^{4}-1=q^{2}>0,$

i.e., parameter $q$ introduced above is positive and finite. Therefore, it
does not depend on parameter $a$. This means that $\mu=1+iq$ is also finite
and does not depend on parameter $a$. Now, taking into account well-known
approximation for the gamma functions

$\frac{{\Gamma\left({z+\alpha}\right)}}{{\Gamma\left({z+\beta}\right)}}\cong
z^{\alpha-\beta},$

one observes that

$\mathop{\lim}\limits_{a\to\infty}\left({\frac{x}{a}+1}\right)^{-\gamma-\lambda_{0}^{2}a^{2}}=\mathop{\lim}\limits_{a\to\infty}\left({\frac{x}{a}+1}\right)^{-\frac{\mu}{2}}=1,$
(4.18)

and

$\mathop{\lim}\limits_{a\to\infty}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{3}}}{{x+a}}}\to\mathrm{e}^{-\lambda_{0}^{2}a^{2}},\quad\mathop{\lim}\limits_{a\to\infty}{\kern
1.0pt}_{1}F_{1}\left({\begin{array}[]{*{20}c}\gamma\\\ \mu\\\
\end{array}\mathchar
24635\relax\;\frac{{2\lambda_{0}^{2}a^{3}}}{{x+a}}}\right)\to{\kern
1.0pt}_{0}F_{1}\left({\begin{array}[]{*{20}c}-\\\ \mu\\\ \end{array}\mathchar
24635\relax\;-2\lambda_{0}^{4}a^{4}}\right).$ (4.19)

Hypergeometric functions ${}_{0}F_{1}$ can be expressed via Bessel functions
of first kind $J_{\alpha}\left(x\right)$ as follows prudnikov2002 :

$_{0}F_{1}\left({\begin{array}[]{*{20}c}-\\\ \mu\\\ \end{array}\mathchar
24635\relax\;-2\lambda_{0}^{4}a^{4}}\right)=\Gamma\left(\mu\right)\left({2\lambda_{0}^{4}a^{4}}\right)^{\frac{{1-\mu}}{{2\lambda_{0}^{4}a^{4}}}}J_{\mu-1}\left({2\sqrt{2}\lambda_{0}^{2}a^{2}}\right).$
(4.20)

The following approximation for Bessel functions of first kind
$J_{\alpha}\left(x\right)$ at $\lvert x\rvert\to\infty$ is well known, too:

$J_{\alpha}\left(x\right)\cong\sqrt{\frac{2}{{\pi x}}}.$ (4.21)

Now, combining together (4.18)-(4.21), one observes that

$\mathop{\lim}\limits_{a\to\infty}\psi_{E}\left(x\right)=0.$

In the end, we would like to discuss briefly the factorization possibilities
of non-relativistic Hamiltonian (2.2) being the sum of BenDaniel–Duke kinetic
energy operator (3.1) and the potential (3.3). In principle, due to that one
replaces constant mass with the mass changing with position, then, the method
for the factorization is different than the method applied for the ordinary
quantum harmonic oscillator, which factorization leads to obtaining its
creation and annihilation operators $\hat{a}^{+}$ and $\hat{a}^{-}$ (2.13).
The method that can be successfully applied here is already known dabrowka1988
; cooper1995 ; plastino1999 ; gonul2002 ; dong2007 ; amir2016 . One writes
down the ground state wavefunction extracted from (3.37) as follows:

$\psi_{0}^{QW}\left(x\right)=C_{0}^{QW}\left({\frac{x}{a}+1}\right)^{-\lambda_{0}^{2}a^{2}}\mathrm{e}^{-\frac{{\lambda_{0}^{2}a^{3}}}{{x+a}}}.$
(4.22)

Next, for simplicity one needs to introduce $\alpha_{0}\left(x\right)$ and
$\rho\left(x\right)$ functions as follows:

$\alpha_{0}\left(x\right)=\frac{{\psi_{0}^{\prime}\left(x\right)}}{{\psi_{0}\left(x\right)}}=-\frac{{\lambda_{0}^{2}a^{2}}}{{x+a}}+\frac{{\lambda_{0}^{2}a^{3}}}{{\left({x+a}\right)^{2}}},\quad\rho\left(x\right)=\frac{{\hbar^{2}}}{{2M}}.$

Then, one can factorize Hamiltonian (2.2) being sum of BenDaniel–Duke kinetic
energy operator (3.1) and the potential (3.3) in terms of operators
$\hat{A}^{+}$ and $\hat{A}^{-}$ as follows:

$\hat{H}=\hbar\omega\left({\hat{A}^{+}\hat{A}^{-}+\frac{1}{2}}\right),$

where,

$\displaystyle\hat{A}^{+}=-\left({\frac{\mathrm{d}}{{\mathrm{d}x}}+\alpha_{0}\left(x\right)}\right)\sqrt{\frac{\rho\left(x\right)}{\hbar\omega}},$
(4.23)
$\displaystyle\hat{A}^{-}=\sqrt{\frac{\rho\left(x\right)}{\hbar\omega}}\left({\frac{\mathrm{d}}{{\mathrm{d}x}}-\alpha_{0}\left(x\right)}\right).$

One can easily prove that the following equation is correct:

$\hat{A}^{-}\psi_{0}^{QW}\left(x\right)=0,$

as well as the following two limit relations also hold:

$\mathop{\lim}\limits_{a\to+\infty}\hat{A}^{\pm}=\hat{a}^{\pm}.$

Only the difference here is that both creation and annihilation operators
$\hat{a}^{\pm}$ together with quantum harmonic oscillator Hamiltonian generate
closed Heisenberg–Weyl algebra and it is dynamical symmetry algebra of the
non-relativistic oscillator model described through Hermite polynomials.
However, here it is impossible to write down the dynamical symmetry algebra of
the model under study by a simpler way just employing $\hat{A}^{\pm}$
operators. Therefore, being focused on the analytical solution of the model
under construction, we leave the problem related to its dynamical symmetry
algebra for our future studies.

We constructed an exactly-solvable model of the non-relativistic harmonic
oscillator with a position-dependent effective mass. It behaves itself as a
semi-infinite quantum well of the non-rectangular profile and its profile is a
step-harmonic potential, but with smooth transition between the step and
harmonic oscillator potentials as a consequence of the certain analytical
dependence of the effective mass from the position and semiconfinement
parameter $a$. We were able to solve the Schrödinger equation corresponding to
the model under construction and obtained that its wavefunctions of the
stationary states are expressed through the ${}_{1}F_{1}$ hypergeometric
functions if the spectrum is continuous and through the Bessel polynomials if
the spectrum is discrete. We also found that the discrete energy spectrum of
the model under construction is non-equidistant and finite, but, at the limit,
when the parameter $a$ goes to infinity, both wavefunctions, and the energy
spectrum correctly recover wavefunctions and the energy spectrum of the usual
non-relativistic harmonic oscillator with a constant effective mass.

We are not aware that a similar exactly-solvable model with the wavefunctions
expressed through the Bessel polynomials has been before constructed.
Therefore, we believe that the exact-solubility advantage of the present model
will make it definitely useful in the description of the various phenomena in
physics.

## Acknowldgement

The authors would like to thank the anonymous reviewer for the valuable
comments and suggestions, which substantially increased the quality of the
paper.

## Conflict of interest

The authors declare that they have no conflicts of interest.

## Data availability

The paper has no associated data, however, any additional information
regarding depicted figures generated from the computations done above is
always available from the corresponding author upon reasonable request.

## References

  * (1) R.C. Miller, A.C. Gossard, D.A. Kleinman and O. Munteanu, Parabolic quantum wells with the $GaAs-Al_{x}Ga_{1-x}As$ system, Phys. Rev. B, 29 3740–-3743 (1984).
  * (2) R.C. Miller, D.A. Kleinman and A.C. Gossard, Energy-gap discontinuities and effective masses for $GaAs-Al_{x}Ga_{1-x}As$ quantum wells, Phys. Rev. B, 29 7085–-7087 (1984).
  * (3) R.C. Miller, A.C. Gossard and D.A. Kleinman, Band offsets from two special $GaAs-Al_{x}Ga_{1-x}As$ quantum well structures, Phys. Rev. B, 32 5443–-5446 (1985).
  * (4) A.C. Gossard, R.C. Miller and W. Wiegmann, MBE growth and energy levels of quantum wells with special shapes, Surf. Sci., 174 131–-135 (1986).
  * (5) L. Rizzi, O.F. Piattella, S.L. Cacciatori and V. Gorini, The step-harmonic potential, Am. J. Phys., 78 842–-850 (2010).
  * (6) A. Amthong, WKB approximation for abruptly varying potential wells, Eur. J. Phys., 35 065009 (2014).
  * (7) J.R. Morris, New scenarios for classical and quantum mechanical systems with position-dependent mass, Quantum Stud.: Math. Found., 2 359–370 (2015).
  * (8) J.R. Morris, Short note: Hamiltonian for a particle with position-dependent mass, Quantum Stud.: Math. Found., 4 295–299 (2017).
  * (9) R. Koekoek, P.A. Lesky and R.F. Swarttouw, Hypergeometric Orthogonal Polynomials and Their $q$-Analogues, Springer, Berlin 2010.
  * (10) P.A.M. Dirac, The quantum theory of the emission and absorption of radiation, Proc. R. Soc. A, 114 243–-265 (1927).
  * (11) L. Infeld and T.E. Hull, The factorization method, Rev. Mod. Phys., 23 21–-68 (1951).
  * (12) A. Messiah, Quantum Mechanics (Vol. I), Wiley, North Holland 1966.
  * (13) D.J. BenDaniel and C.B. Duke, Space-charge effects on electron tunneling, Phys. Rev., 152 683–-692 (1966).
  * (14) W.A. Harrison, Tunneling from an independent-particle point of view, Phys. Rev., 123 85–-89 (1961).
  * (15) I. Giaever, Energy gap in superconductors measured by electron tunneling, Phys. Rev. Lett., 5 147–-148 (1960).
  * (16) I. Giaever, Electron tunneling between two superconductors, Phys. Rev. Lett., 5 464–-466 (1960).
  * (17) T. Gora and F. Williams, Theory of electronic states and transport in graded mixed semiconductors, Phys. Rev., 177 1179–-1182 (1969).
  * (18) Q.-G. Zhu and H. Kroemer, Interface connection rules for effective-mass wave functions at an abrupt heterojunction between two different semiconductors, Phys. Rev. B, 27 3519–-3527 (1983).
  * (19) O. von Roos, Position-dependent effective masses in semiconductor theory, Phys. Rev. B, 27 7547–-7551 (1983).
  * (20) R.A. Morrow and K.R. Brownstein, Model effective-mass Hamiltonians for abrupt heterojunctions and the associated wave-function-matching conditions, Phys. Rev. B, 30 678–-680 (1984).
  * (21) T.L. Li and K.J. Kuhn, Band-offset ratio dependence on the effective-mass Hamiltonian based on a modified profile of the $GaAs$ \- $Al_{x}Ga_{1-x}As$ quantum well, Phys. Rev. B, 47 12760–-12770 (1993).
  * (22) J. Li and M. Ostoja-Starzewski, Fractal solids, product measures and fractional wave equations, Proc. R. Soc. A, 465 2521–-2536 (2009).
  * (23) J.R.F. Lima, M. Vieira, C. Furtado, F. Moraes and C. Filgueiras, Yet another position-dependent mass quantum model, J. Math. Phys., 53 072101 (2012).
  * (24) F.D. Nobre, M.A. Rego-Monteiro, Non-Hermitian PT Symmetric Hamiltonian with position-dependent masses: associated Schrödinger equation and finite-norm solutions, Braz. J. Phys., 45, 79–88 (2015).
  * (25) O. Mustafa, Position-dependent mass momentum operator and minimal coupling: point canonical transformation and isospectrality, Eur. Phys. J. Plus, 134, 228 (2019).
  * (26) R.A. El-Nabulsi, On a new fractional uncertainty relation and its implications in quantum mechanics and molecular physics, Proc. R. Soc. A, 476 20190729 (2020).
  * (27) A.V. Kolesnikov and A.P. Silin, Quantum mechanics with coordinate-dependent mass, Phys. Rev. B, 59, 7596–7599 (1999).
  * (28) A.F. Nikiforov, and V.B. Uvarov, Special Functions of Mathematical Physics: A Unified Introduction with Applications, Birkhäuser, Basel 1988.
  * (29) L.D. Landau and E.M. Lifshitz, Quantum Mechanics: Non-Relativistic Theory, Pergamon Press, Oxford 1991.
  * (30) E.I. Jafarov and J. Van der Jeugt, Exact solution of the semiconfined harmonic oscillator model with a position-dependent effective mass, Eur. Phys. J. Plus, 136 758 (2021).
  * (31) E.I. Jafarov and J. Van der Jeugt, Exact solution of the semiconfined harmonic oscillator model with a position-dependent effective mass in an external homogeneous field, Pramana – J. Phys., 96, 35 (2022).
  * (32) C. Quesne, Generalized semiconfined harmonic oscillator model with a position-dependent effective mass, Eur. Phys. J. Plus, 137, 225 (2022).
  * (33) P.A. Lesky, Einordnung der Polynome von Romanovski-Bessel in das Askey-Tableau, Z. Angew. Math. Mech., 78, 646–648 (1998).
  * (34) E.J. Routh, On some properties of certain solutions of a differential equation of the second order, Proc. Lond. Math. Soc., 16, 245–261 (1885).
  * (35) V.I. Romanovski, Sur quelques classes nouvelles de polynomes orthogonaux, C. R. Acad. Sci. Paris, 188, 1023–1025 (1929).
  * (36) E.I. Jafarov, A.M. Mammadova and J. Van der Jeugt, On the direct limit from pseudo Jacobi polynomials to Hermite polynomials, Mathematics, 9, 88 (2021).
  * (37) S.M. Nagiyev, On two direct limits relating pseudo-Jacobi polynomials to Hermite polynomials and the pseudo-Jacobi oscillator in a homogeneous gravitational field, Theor. Math. Phys., 210, 121–134 (2022).
  * (38) A.P. Prudnikov, Y.A. Brychkov and O.I. Marichev, Integrals and Series: vol.3 – More Special Functions, Taylor and Francis, London 2002.
  * (39) J.W. Dabrowska, A. Khare and U.P. Sukhatme, Explicit wavefunctions for shape-invariant potentials by operator techniques, J. Phys. A: Math. Gen., 21, L195–L200 (1988).
  * (40) F. Cooper, A. Khare and U. Sukhatme, Supersymmetry and quantum mechanics, Phys. Rep., 251, 267–385 (1995).
  * (41) A.R. Plastino, A. Rigo, M. Casas, F. Garcias and A. Plastino, Supersymmetric approach to quantum systems with position-dependent effective mass, Phys. Rev. A, 60, 4318–4325 (1999).
  * (42) B. Gönül, B. Gönül, D. Tutcu and O. Özer, Supersymmetric approach to exactly solvable systems with position-dependent effective masses, Mod. Phys. Lett. A, 17, 2057–2066 (2002).
  * (43) S.-H. Dong, J.J. Peña, C. Pachego-García and J. García-Ravelo, Algebraic approach to the position-dependent mass Schrödinger equation for a singular oscillator, Mod. Phys. Lett. A, 22, 1039–1045 (2007).
  * (44) N. Amir and S. Iqbal, Algebraic solutions of shape-invariant position-dependent effective mass systems, J. Math. Phys., 57, 062105 (2016).


Long<EMAIL_ADDRESS>Yuli<EMAIL_ADDRESS>Dorit Merhofdorit.merhof@ur.de2,3 Institute of Imaging & Computer Vision
RWTH Aachen University
Aachen, Germany Faculty of Informatics and Data Science, University of
Regensburg, Germany Fraunhofer Institute for Digital Medicine MEVIS, Bremen,
Germany Instance Segmention via Layering

# Instance Segmentation of Dense and Overlapping Objects via Layering

###### Abstract

Instance segmentation aims to delineate each individual object of interest in
an image. State-of-the-art approaches achieve this goal by either partitioning
semantic segmentations or refining coarse representations of detected objects.
In this work, we propose a novel approach to solve the problem via object
layering, i.e. by distributing crowded, even overlapping objects into
different layers. By grouping spatially separated objects in the same layer,
instances can be effortlessly isolated by extracting connected components in
each layer. In comparison to previous methods, our approach is not affected by
complex object shapes or object overlaps. With minimal post-processing, our
method yields very competitive results on a diverse line of datasets: C.
elegans (BBBC), Overlapping Cervical Cells (OCC) and cultured neuroblastoma
cells (CCDB). The source code is publicly available
111https://github.com/looooongChen/instSeg/.

##  1 Introduction

Different from semantic segmentation, which pays no attention to the
individual objects, instance segmentation aims not only to associate every
pixel of an image with a class label but also to delineate objects of the same
class as individuals. This task becomes more challenging, when objects are
densely located or even overlapping with each other.

A prevalent top-down class of instance segmentation methods is detection-
based, which firstly localize an instance with a coarse shape representation
and then refine the shape in an additional step. As a paradigm, Mask-RCNN [He
et al.(2018)He, Gkioxari, Dollár, and Girshick] refines bounding boxes
obtained from Region-based Convolutional Neural Networks (R-CNN) [Girshick et
al.(2014)Girshick, Donahue, Darrell, and Malik, Ren et al.(2015)Ren, He,
Girshick, and Sun]. Relying on non-maximum suppression (NMS) to remove
duplicate predictions, detection-based methods become less competent in the
cases of dense clusters and overlapping objects. A finer polygon
representation approach was proposed by [Schmidt et al.(2018)Schmidt, Weigert,
Broaddus, and Myers] to reduce false suppressions. However, NMS has a
methodological flaw for objects that inherently overlap. Moreover, many
objects, especially in the biomedical domain, cannot be well-approximated by
tractable shape representations, such as bounding boxes and star-convex
polygons.

Free from false suppression and coarse shape approximation, alternative
bottom-up approaches obtain instances by grouping pixels [Chen et
al.(2016)Chen, Qi, Yu, and Heng, Gao et al.(2019)Gao, Shan, Wang, Zhao, Yu,
Yang, and Huang, Bai and Urtasun(2017)]. DCAN [Chen et al.(2016)Chen, Qi, Yu,
and Heng] predicts the object boundary explicitly and groups connected pixels
that are separated by these boundaries as instances. However, this approach is
sensitive to broken boundaries. A few misclassified pixels can lead to
erroneous merging of adjacent objects. Graph partition based on pixel-pair
affinity [Gao et al.(2019)Gao, Shan, Wang, Zhao, Yu, Yang, and Huang] and the
watershed transform [Bai and Urtasun(2017)] are more robust in terms of
grouping, but the pipeline as a whole relies heavily on post-processing, with
the learning model only optimized with intermediate results. In addition,
grouping-based approaches are inherently incapable of handling overlapping
objects, where one pixel may belong to more than one object.

Recent research [De Brabandere et al.(2017)De Brabandere, Neven, and Van Gool,
Chen et al.(2019)Chen, Strauch, and Merhof, Kong and Fowlkes(2018), Wu et
al.(2020)Wu, Chen, and Merhof] introduces pixel embedding for grouping. In
these approaches, a deep neural network is trained to map pixels into an
embedding space, in which pixel embeddings from the same object are close,
while those from different objects, especially adjacent objects, are apart.
Then, pixels are grouped in the embedding space with low-level clustering
algorithms, such as Mean Shift [Comaniciu and Meer(2002)] and DBSCAN [Ester et
al.(1996)Ester, Kriegel, Sander, Xu, et al.]. Although pixel embedding based
grouping has proven to be more robust, it still can not avoid the shortcomings
of grouping-based methods: the reliance on post-processing, the optimization
of intermediate results, and the incapability to segment overlapping objects.

Our work is partially inspired by the pixel embedding approach. To alleviate
the aforementioned limitations, we propose to train a more structured
embedding space. Specifically, we restrict one object to "live" in one
dimension. Correspondingly, embedding vectors will be one-hot in overlap-free
areas and have more than one active digits at locations where an overlap of
objects occurs (Figure 1). We figuratively call this process of objects being
distributed to different dimensions "layering", and use the terms "dimension"
and "layer" interchangeably in the following context. Since our training loss
penalizes adjacent objects having the same embedding vector, only spatially
separated objects can be assigned to the same layer, eliminate crowding and
overlapping of objects. Therefore, instances can be effortlessly obtained by
finding connected components in each layer.

The main contributions of our work are as follows: We propose (1) a novel
approach to eliminate crowdedness and overlap of instances by layering objects
into different output layers; (2) an approach of _spontaneous_ object layering
through deep model learning; (3) a concise and effective framework for
segmentation of densely distributed objects without data-specific post-
processing efforts.

To our best knowledge, our work is the first pixel embedding based approach
that does not require explicit pixel clustering, and is capable of handling
object overlap. Competing with several state-of-the-art approaches, our method
yields comparable or better results on a diverse line of datasets: C. elegans
[Ljosa et al.(2012)Ljosa, Sokolnicki, and Carpenter] (BBBC), overlapping
cervical cells [Lu et al.(2015)Lu, Carneiro, and Bradley, Lu et al.(2016)Lu,
Carneiro, Bradley, Ushizima, Nosrati, Bianchi, Carneiro, and Hamarneh] (OCC)
and cultured neuroblastoma cells [Yu et al.()Yu, Lee, Hariharan, Bu, and
Ahmed] (CCDB).

Figure 1: Our approach eliminates the object crowding and overlapping by
distributing adjacent objects into different layers. (a) To achieve this, each
foreground pixel is assigned by drawing pixels of the same object into the
same layer and pushing pixels of adjacent objects into different layers
(Section 2.2). Areas where overlap occurs are ignored in this training phase.
(b) After object layering converges, overlapping areas are trained to complete
the object (Section 2.3). At locations where overlap occurs (the red cross and
line), more than one layer will give foreground prediction.

##  2 Proposed Method

###  2.1 Overview

Our model consists of two output branches: the foreground branch and the
layering branch (Figure 1). The foreground segmentation, as an auxiliary,
excludes background pixels from further processing, while the layering branch
is devoted to separating foreground instances.

The foreground branch uses a 1-channel convolutional layer with _Sigmoid_
activation for output. The output layer can also be set to multi-channel and
_Softmax_ activated in the case of more than one semantic category. The
layering branch has $N$ output channels with the _Sigmoid_ activation, each of
which is trained to only contain spatially separated objects. Since objects in
the same layer exhibit no contact or overlap, they can be effortlessly
isolated by extracting connected components.

We train the model in two phases: layering training and overlap completion. In
the first phase, under the supervision of our proposed layering loss (Section
2.2), the model learns to assign foreground pixels to one of the layers,
maintaining the restriction that neighboring objects should be located in
different layers. We only consider image areas where no objects overlap in
this phase, since the layering loss exclusively chooses one layer for one
pixel. After the layering training converges, we generate $N$ binary masks by
ordering each object mask into one of the $N$ layers, based on the layering
results of the object’s overlap-free part. It is worth mentioning that more
than one layer will be positive at locations where an overlap of objects
occurs. A _Dice_ -like loss (Section 2.3) is then computed with the generated
masks and further included in the second training phase. The overlapping area
is explicitly trained in this phase, to complete the intact object.

In summary, the model is trained in two phases with the following loss:

$L=[L_{foreground}]_{S}^{1,2}+[L_{layering}]_{S_{fn}}^{1,2}+[L_{overlap}]_{S_{f}}^{2},$
(1)

where the subscript of $[\cdot]$ denotes on which area a loss term is
computed: $S$, $S_{f}$ and $S_{fn}$ represent three progressively smaller
areas, namely, the whole image, the foreground and the foreground without
object overlap. The superscript indicates in which training phase the loss
term is included. The standard _Crossentropy_ is used as foreground training
loss $L_{foreground}$. Details of $L_{layering}$ and $L_{overlap}$ are
depicted in the following two sections.

###  2.2 Layering Loss

The object layering is achieved by drawing pixels of the same object together
into the same dimension with an attracting loss term $L_{attr}$ and pushing
away neighboring objects into different dimensions with a repelling loss term
$L_{rep}$. Our loss is constructed aroud the cosine similarity
$D(\mathbf{e_{i}},\mathbf{e_{j}})=\frac{\mathbf{e_{i}}^{\top}\mathbf{e_{j}}}{\|\mathbf{e_{i}}\|_{2}\|\mathbf{e_{j}}\|_{2}}$
[Chen et al.(2019)Chen, Strauch, and Merhof, Wu et al.(2020)Wu, Chen, and
Merhof], which is zero when two vectors are located in two othogonal spaces.
The operation $\|\cdot\|_{2}$ computes the $\mathcal{L}^{2}$ norm.

We push adjacent objects into the orthogonal space of each other, while non-
adjacent objects can stay in the same cluster. Assuming that there are $C$
objects ($\\{O_{i}|i=1,2,...,C\\}$) in an image, we represent the overlapping
and overlap-free part of an object with $O_{i}^{o}$ and $O_{i}^{n}$,
respectively. The attracting and repelling terms can thus be formulated as:

$\displaystyle L_{attr}=$ $\displaystyle
1-\frac{1}{\sum_{i=1}^{C}|O_{i}^{n}|}\sum_{i=1}^{C}\sum_{p\in
O_{i}^{n}}D(\mathbf{e_{p}},\mathbf{u_{i}})^{2}\;,$ (2) $\displaystyle
L_{rep}=$
$\displaystyle\frac{1}{C}\sum_{i=1}^{C}\frac{1}{|Adj(O_{i})|}\sum_{j\in
Adj(O_{i})}D(\mathbf{u_{i}},\mathbf{u_{j}})^{2}\;,$ (3)

where $e_{p}$ indicates the embedding vector of pixel $p$ and
$\mathbf{u_{i}}=\frac{1}{|O_{i}^{n}|}\sum_{p\in O_{i}^{n}}\mathbf{e_{p}}$ is
the mean embedding of overlap-free part of the $i$-th object. $Adj(O_{i})$
represents the set of adjacent objects to object $O_{i}$, whose shortest
distances to $O_{i}$ are less than a threshold $t$ ($t=15$ pixels in this
work). The operator $|\cdot|$ returns the element number of a set. For
example, $|O_{i}^{n}|$ is the pixel number of area $O_{i}^{n}$ and
$|Adj(O_{i})|$ is the number of adjacencies of object $O_{i}$.

Training with the attracting term $L_{attr}$ and the repelling term $L_{rep}$
specifies embedding vectors to locate in mutually orthogonal spaces, but it
does not guarantee that the vectors only "live" in one layer (dimension), i.e.
they are distributed on the standard axes (see Figure 2). Therefore, we
introduce the sparse term $L_{sparse}$, which imposes preferences for the
vector whose digits, except for a single one, are suppressed to 0:

$\displaystyle L_{sparse}=$ $\displaystyle
1-\frac{1}{\sum_{i=1}^{C}|O_{i}^{n}|}\sum_{i=1}^{C}\sum_{p\in
O_{i}^{n}}\max(\frac{\mathbf{e_{p}}}{||\mathbf{e_{p}}||_{2}})\;,$ (4)

where the operator $\max(\cdot)$ takes the value of the maximal digit.

The layering loss consists of the three terms above:

$\displaystyle L_{layering}=L_{attr}+L_{rep}+\lambda L_{sparse}\;,$ (5)

where $\lambda$ is a weighting constant, and we use $\lambda=0.1$ in this
work. It is worth mentioning that only the overlap-free foreground part is
involved in the calculation of $L_{layering}$.

Figure 2: Visualization of trained embedding space. (a) Four objects with
three overlapping areas are present in this example. (b) After layering
training, the object parts that do not overlap are ideally distributed in
different dimensions (layers). Spatially separated objects, such as A and D,
can exist in the same dimension. Pixels of overlapping areas, untrained in
this phase, live randomly in the plane composed of axes to which the
overlapping objects belong. (c) After overlap completion, pixels of the
overlapping area will congregate to the point $(1,1)$, indicating a foreground
area in both layers. And, pixels of the non-overlap part are also more
compactly distributed and closer to 1.

###  2.3 Overlap Completion

After layering training, the overlap-free parts of objects are assigned to one
of the layers, while the overlapping parts remain untrained. Generally, the
untrained pixels from overlapping areas will randomly belong to one dimension
to which one of the objects that overlap belongs (Figure 2b). Accordingly, a
random boundary between objects will be formed, as illustrated in Figure 1a.

To train the region with overlap, we firstly generate $N$ binary masks,
denoted by $S\in\\{0,1\\}^{H\times W\times N}$, by placing the silhouette of
each object into one of the $N$ masks based on the layering results of the
object’s overlap-free part. Indexing the stack of N layers at different pixel
positions $p$ with $S_{p}\in\\{0,1\\}^{N}$, we define the following _Dice_
-like [Milletari et al.(2016)Milletari, Navab, and Ahmadi] loss for overlap
completion training:

$\displaystyle L_{overlap}=$ $\displaystyle 1-2\frac{\sum_{i=1}^{C}\sum_{p\in
O_{i}}\mathbf{e_{p}}^{\top}\mathbf{S_{p}}}{\sum_{i=1}^{C}\sum_{p\in
O_{i}}(\mathbf{1}^{\top}\mathbf{e_{p}}+\mathbf{1}^{\top}\mathbf{S_{p}})}\;.$
(6)

The dot product between binary vector $\mathbf{e_{p}}$ and $\mathbf{S_{p}}$ on
the numerator of $L_{overlap}$ is analogous to the intersection in _Dice_
loss, whose value is maximized when $\mathbf{e_{p}}$ and $\mathbf{S_{p}}$ are
exactly the same. The denominator sums up all digits of vector
$\mathbf{e_{p}}$ and $\mathbf{S_{p}}$, which is equivalent to the sum in the
_Dice_ loss. The overlap completion loss $L_{overlap}$ is applied to the
foreground, including the overlapping and non-overlapping areas. The vector
$\mathbf{S_{p}}$ is one-hot in areas without overlap, while it has more than
one non-zero digits at locations where overlap occurs.

###  2.4 Post-Processing

Since each of the $N$ layers only contains spatially separated objects, high-
quality instance segmentations can be obtained with minimal post-processing
effort, only involving simple and computationally efficient operations. Our
post-processing steps are also general for all instance segmentation tasks,
without utilizing any prior knowledge specific to a certain dataset. Detailed
steps are listed in Algorithm 1. It worth mentioning that line 7-8 can be
ignored if there is no object overlap for certain tasks.

The post-processing requires two trivial parameters: a threshold $\tau$ and a
minimal object size $S_{min}$. The value $\tau$ is used for foreground
thresholding (lines 1 and 7 in Algorithm 1). The minimal object size $S_{min}$
is responsible for eliminating small noisy objects (lines 2 and 11 in
Algorithm 1). We use $\tau=0.5$ and $S_{min}=250$ for all experiments in this
work.

Algorithm 1 Post-processing to obtain object segmentations

1:(raw prediction) foreground $F_{raw}\in(0,1)^{H\times W}$, layering
$L_{raw}\in(0,1)^{H\times W\times N}$

2:(parameters) threshold $\tau$, minimal object size $S_{min}$

3:threshold $F_{raw}$ with value $\tau$ to get foreground $F$

4:remove connected components in $F$ that are smaller than $S_{min}$

5:initial $N$ empty binary layering masks $L\in\\{0,1\\}^{H\times W\times N}$

6:for each pixel location $i\in\\{1,2,...,H\\},j\in\\{1,2,...W\\}$ and layer
$k\in\\{1,2,...,N\\}$ do

7: if $L_{raw}(i,j,k)$ is the largest in $L_{raw}(i,j,:)$ then

8: set $L(i,j,k)$ to 1

9: else if $L_{raw}(i,j,k)$ is larger than $\tau$ then $\triangleright$ omit
if no overlap exists

10: set $L(i,j,k)$ to 1 $\triangleright$ omit if no overlap exists

11: end if

12:end for

13:take all connected components larger than $S_{min}$ in each layer as
objects

##  3 Experiments and Results

###  3.1 Implementation

While our approach is not tied to a particular network architecture, we
perform experiments with UNet [Ronneberger et al.(2015)Ronneberger, Fischer,
and Brox]. In our implementation, we apply batch normalization after each
convolutional layer. In addition, we experiment with a UNet-S and a UNet-L,
which have, respectively, one less and one more convolutional block on the
contracting and expansive path than the original UNet (denoted as UNet-M), to
investigate the effect of receptive field size.

The foreground branch and the layering branch share the same feature map from
the UNet backbone. Two 3x3 covolutional layers with 64 features are added
before the output layer to avoid "feature conflict". We use an $N$-channel
convolutional layer ($N=8$ in this work) and a $1$-channel convolutional layer
with the _sigmoid_ activation for the final output.

The model was trained with _RMSprop_ optimizer [Hinton et al.(2012)Hinton,
Srivastava, and Swersky] with a learning rate of $1e^{-4}$ exponentially
decayed to 0.9 every 10000 steps. The training dataset was randomly split into
90$\%$ and 10$\%$ for training and validation. The layering training phase
lasted 1000 epoches. Then, the overlap completion training continued from the
best model of layering training for another 500 epoches. The best model of
overlap completion training was used as the final model for evaluation. The
"best" models are chosen based on validation results.

###  3.2 Datasets and Evaluation Metrics

For the evaluation, we use three biomedical image datasets containing a
population of objects with different shapes, different degrees of density and
overlap. All datasets were augmented using random horizontal and vertical
flips, random rotation, random gamma $\gamma\in(0.5,2)$ correction transform
and elastic deformation.

BBBC010 contains 100 bright-field microscopic images of live/dead C. elegans
[Ljosa et al.(2012)Ljosa, Sokolnicki, and Carpenter]. The C. elegans are
slender, bilaterally symmetrical objects in curved or ring-shaped poses, which
may cross over others in this data set. Subset D was left for evaluation and
the rest for training. We crop the image and a 448x448 pixel area containing
objects remains.

OCC2014 is an EDF (extended depth of field) image collection of overlapping
cervical cells from Pap smears [Lu et al.(2015)Lu, Carneiro, and Bradley, Lu
et al.(2016)Lu, Carneiro, Bradley, Ushizima, Nosrati, Bianchi, Carneiro, and
Hamarneh], consisting of 16 real and 945 synthesized images. The cells are
roundish and densely clustered, and the overlap can be larger than 50$\%$ of
the object size. In our experiments, the original training set and half of the
test set served for training. Evaluation was conducted on the other half. The
images are resized to 320x320 pixels.

CCDB6843 contains 100 wide field fluorescent images of cultured neuroblastoma
cells collected by [Yu et al.()Yu, Lee, Hariharan, Bu, and Ahmed]. The cells
are densely located and irregularly shaped. We randomly chose 24 images as the
test set. All images are resized to 448x448 pixels.

In our evaluation, a predicted object $I_{pred}$ is considered to be a
successful match (true positive $TP_{t}$) if its intersection over union
$IoU=\frac{I_{pred}\cap I_{gt}}{I_{pred}\cup I_{gt}}$ with a ground truth
object is larger than a given threshold $t$, while unmatched predictions and
ground truth objects are counted as false positive ($FP_{t}$) and false
negative ($FN_{t}$). Using these values, a measure of detection accuracy can
be calculated: $AP_{t}=\frac{TP_{t}}{TP_{t}+FP_{t}+FN_{t}}$. By passing from
loose to strict thresholds, the segmentation accuracy is also reflected. To
better quantify the segmentation performance, we also calculate the Aggregated
Jaccard Index ($AJI$) [Kumar et al.(2017)Kumar, Verma, Sharma, Bhargava,
Vahadane, and Sethi].

Layer-1

Layer-2

Layer-3

Layer-4

Layer-5

Layer-6

Layer-7

Layer-8

Figure 3: Activation map of each layer: results on the datasets BBBC, OCC and
CCDB from top to bottom. Spatially adjacent objects are distributed to
different layers, eliminating crowding and overlapping of objects. Meanwhile,
not all layers contain objects, which suggests that, by keeping only adjacent
objects separate, layers are fully utilized.

Image

UNet-2

UNet-3

M-RCNN

StarDist

Ours

GT

Figure 4: Qualitative segmentation results on the datasets BBBC, OCC and CCDB
from top to bottom. A few typical errors are marked with arrows: false
suppression (red), merged objects (green), inaccurate shape (yellow) and
incapability to handle overlap (blue).

###  3.3 Competing Methods

UNet-2/3: We train UNet [Ronneberger et al.(2015)Ronneberger, Fischer, and
Brox] models as semantic segmentation tasks with a two labels (object and
background) and a three labels (object, boundary and background) setup. In
terms of network structure, UNet is the least different from our model: UNet
output all objects with a single image plane, whereas our model has several
containing layered objects.

Mask-RCNN: Mask-RCNN [He et al.(2018)He, Gkioxari, Dollár, and Girshick]
localizes objects by proposal classification and non-max suppression (NMS).
Then, segmentation is performed on each object bounding box. The NMS threshold
was set to 0.7 for all experiments.

StarDist: Without an explicit segmentation step, StarDist [Schmidt et
al.(2018)Schmidt, Weigert, Broaddus, and Myers] obtains object masks by
combining distances from the center to the boundary along different
directions. We used 32 radial directions in our experiments.

All methods except Mask-RCNN were trained from scratch, while the Mask-RCNN
model was fine-tuned on a COCO pretrained model [Lin et al.(2014)Lin, Maire,
Belongie, Hays, Perona, Ramanan, Dollár, and Zitnick] with the ResNet-101
backbone [He et al.(2016)He, Zhang, Ren, and Sun]. For the 3-label UNet and
StarDist, we used a pseudo-boundary as an approximated separation of
overlapping regions, obtained by skeletonizing [Lee et al.(1994)Lee, Kashyap,
and Chu] the binary mask of overlapping regions together with boundaries.

###  3.4 Results and Discussion

For a careful interpretation, we discuss the methods from the following
aspects: the handling of object shape, touching objects, object overlap and
the effect of high object density.

Shape: Using radial directions, StarDist has difficulty in reconstructing
slender shapes and boundaries with fine structure. This is reflected by the
broken body of C. elegans and the roundish approximation of the neuroblastoma
cells (Figure 4). Other methods, which perform pixel-wise segmentation, are
less affected.

Touching objects: Since it does not take clustered objects explicitly into
consideration, the 2-label UNet performs suboptimally on all datasets due to
false fusion, while others are object-aware to different degrees under
different conditions.

Object overlap: Since UNet and StarDist naturally assign one pixel to one
object, an overlapping area can only be handled with approximated boundaries.
In case of less severe overlap, such as on the BBBC and CCDB datasets, UNet
(3-label) still achieves very good results, although its performances drops
significantly when the overlap ratio increases on the dataset OCC.

Object density: Mask-RCNN suffers from the false suppression of NMS, when the
bounding boxes of two objects overlap with a large ratio, as, for example, in
the case of the two parallel and close C. elegans in Figure 4. We also find
that Mask-RCNN has difficulty in distinguishing two closely located irregular-
shaped objects, such as the neuroblastoma cells in Figure 4.

As described above, segmentation errors, including false suppression, merged
objects, inaccurate shape and the incapability to handle overlap (Figure 4),
occur on different approaches, depending on the data characteristics. By
contrast, our method is more robust in all of these aspects, and, therefore,
achieves the best and, evidently, better results on the dataset BBBC and CCDB
(Table 1). On the dataset OCC, the performance of our approach is only
marginally worse than Mask-RCNN in terms of the Aggregated Jaccard Index
($AJI$).

Table 1: Quantitative evaluation. Average precisions ($AP_{t}$) under different $IoU$ thresholds, mean values of average precision $AP$ over these thresholds, and the Aggregated Jaccard Index ($AJI$) are reported. The best two results are shown in bold, and the best is underlined. Data and Methods | $AP_{0.5}$ | $AP_{0.6}$ | $AP_{0.7}$ | $AP_{0.8}$ | $AP_{0.9}$ | $meanAP$ | $AJI$
---|---|---|---|---|---|---|---
BBBC | UNet-2 | .5455 | .4645 | .4497 | .4212 | .2355 | .4233 | .5327
UNet-3 | .8863 | .8197 | .7412 | .5795 | .2717 | .6597 | .7786
StarDist | .3098 | .1410 | .0372 | .0011 | .000 | .0978 | .4499
MRCNN | .8953 | .8629 | .8111 | .5305 | .0382 | .6276 | .7580
Ours | .9357 | .9188 | .8648 | .7606 | .2904 | .7541 | .8442
OCC | UNet-2 | .1548 | .1303 | .1129 | .1055 | .1048 | .1217 | .2058
UNet-3 | .7010 | .6071 | .5097 | .3767 | .1950 | .4779 | .5802
StarDist | .6556 | .5566 | .4346 | .2970 | .1547 | .4197 | .6927
MRCNN | .9277 | .9181 | .8870 | .8117 | .5564 | .8202 | .8412
Ours | .9230 | .8768 | .8007 | .6788 | .4349 | .7429 | .8353
CCDB | UNet-2 | .3698 | .3360 | .3049 | .2763 | .2228 | .3020 | .1185
UNet-3 | .7307 | .6774 | .6210 | .5153 | .2838 | .5656 | .7148
StarDist | .7428 | .6532 | .4958 | .2685 | .0326 | .4386 | .6903
MRCNN | .6248 | .5691 | 4888 | 3476 | .0763 | .4213 | .5842
Ours | .7968 | .7467 | .6767 | .4889 | .2230 | .5864 | .7601

By chance, through some examples from experiments, we found certain layers
reveal obvious semantics. For example, common morphological features, such as
body orientation, are observed in layers predicted by the BBBC models (Figure
5). Another feature commonly exploited by all models is the relative position.
For example, the leftmost protruding objects are active in one common layer of
OCC2014 predictions (Figure 5). To capture high-level morphological features
and sophisticated relative positions, especially in very crowded cases, an
adequately large receptive field (RF) is required. To verify our analysis, we
train three UNet variants with RF of 108, 220 and 444 pixels (Section 3.1).
The experiments shows very significant differences: 39.49$\%$, 54.29$\%$ and
12.04$\%$ improvement from UNet-S to UNet-L on the datasets BBBC, OCC and CCDB
in terms of $meanAP$ (Table 2).

In addition, we compared the performance before and after overlap completion
(Table 2). On the dataset BBBC and CCDB, the difference is relatively small:
5.91$\%$ and 2.11$\%$, taking the best model UNet-L as an example, since the
objects are only slightly overlapped. By contrast, the segmentation of
severely overlapping cervical cells gains a 36.54$\%$ performance boost.

Table 2: Performance of models with different sizes with/without overlap completion (OC)  $meanAP$ | UNet-S | UNet-M | UNet-L
---|---|---|---
BBBC | w/o OC | .4655 | .6670 | .7120
w/ OC | .5406 | .7162 | .7541
OCC | w/o OC | .4359 | .5273 | .5441
w/ OC | .4815 | .6841 | .7429
CCDB | w/o OC | .5229 | .5569 | .5743
w/ OC | .5234 | .5614 | .5864

Figure 5: Example activation maps of certain layers. _C. elegans_ with the
vertical body orientation are active in the 4th layer of the BBBC data (first
row). In the 8th layer of OCC data (second row), most leftmost protruding
cells are active.

##  4 Conclusion and Outlook

Our proposed approach can successfully layer touching and overlapping objects
into different image layers. By grouping spatially separated objects in the
same layer, our method simplifies post-processing and improves the accuracy of
instance segmentation, yielding very competitive results on a diverse line of
data sets. Our future work will focus on understanding the layering mechanism.

## References

  * [Bai and Urtasun(2017)] Min Bai and Raquel Urtasun.  Deep watershed transform for instance segmentation.  In _Proceedings of the IEEE conference on computer vision and pattern recognition_ , pages 5221–5229, 2017.
  * [Chen et al.(2016)Chen, Qi, Yu, and Heng] Hao Chen, Xiaojuan Qi, Lequan Yu, and Pheng-Ann Heng.  Dcan: deep contour-aware networks for accurate gland segmentation.  In _Proceedings of the IEEE conference on Computer Vision and Pattern Recognition_ , pages 2487–2496, 2016.
  * [Chen et al.(2019)Chen, Strauch, and Merhof] Long Chen, Martin Strauch, and Dorit Merhof.  Instance segmentation of biomedical images with an object-aware embedding learned with local constraints.  In _International Conference on Medical Image Computing and Computer-Assisted Intervention_ , pages 451–459. Springer, 2019.
  * [Comaniciu and Meer(2002)] Dorin Comaniciu and Peter Meer.  Mean shift: A robust approach toward feature space analysis.  _IEEE Transactions on pattern analysis and machine intelligence_ , 24(5):603–619, 2002.
  * [De Brabandere et al.(2017)De Brabandere, Neven, and Van Gool] Bert De Brabandere, Davy Neven, and Luc Van Gool.  Semantic instance segmentation for autonomous driving.  In _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition Workshops_ , pages 7–9, 2017.
  * [Ester et al.(1996)Ester, Kriegel, Sander, Xu, et al.] Martin Ester, Hans-Peter Kriegel, Jörg Sander, Xiaowei Xu, et al.  A density-based algorithm for discovering clusters in large spatial databases with noise.  In _Kdd_ , volume 96, pages 226–231, 1996.
  * [Gao et al.(2019)Gao, Shan, Wang, Zhao, Yu, Yang, and Huang] Naiyu Gao, Yanhu Shan, Yupei Wang, Xin Zhao, Yinan Yu, Ming Yang, and Kaiqi Huang.  Ssap: Single-shot instance segmentation with affinity pyramid.  In _Proceedings of the IEEE/CVF International Conference on Computer Vision_ , pages 642–651, 2019.
  * [Girshick et al.(2014)Girshick, Donahue, Darrell, and Malik] Ross Girshick, Jeff Donahue, Trevor Darrell, and Jitendra Malik.  Rich feature hierarchies for accurate object detection and semantic segmentation, 2014.
  * [He et al.(2016)He, Zhang, Ren, and Sun] Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In _Proceedings of the IEEE conference on computer vision and pattern recognition_ , pages 770–778, 2016.
  * [He et al.(2018)He, Gkioxari, Dollár, and Girshick] Kaiming He, Georgia Gkioxari, Piotr Dollár, and Ross Girshick.  Mask r-cnn, 2018.
  * [Hinton et al.(2012)Hinton, Srivastava, and Swersky] Geoffrey Hinton, Nitish Srivastava, and Kevin Swersky.  Neural networks for machine learning lecture 6a overview of mini-batch gradient descent.  _Cited on_ , 14(8):2, 2012.
  * [Kong and Fowlkes(2018)] Shu Kong and Charless C Fowlkes.  Recurrent pixel embedding for instance grouping.  In _Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition_ , pages 9018–9028, 2018.
  * [Kumar et al.(2017)Kumar, Verma, Sharma, Bhargava, Vahadane, and Sethi] Neeraj Kumar, Ruchika Verma, Sanuj Sharma, Surabhi Bhargava, Abhishek Vahadane, and Amit Sethi.  A dataset and a technique for generalized nuclear segmentation for computational pathology.  _IEEE transactions on medical imaging_ , 36(7):1550–1560, 2017.
  * [Lee et al.(1994)Lee, Kashyap, and Chu] Ta-Chih Lee, Rangasami L. Kashyap, and Chong-Nam Chu.  Building skeleton models via 3-d medial surface/axis thinning algorithms.  _Computer Vision, Graphics, and Image Processing_ , 56(6):462–478, 1994.
  * [Lin et al.(2014)Lin, Maire, Belongie, Hays, Perona, Ramanan, Dollár, and Zitnick] Tsung-Yi Lin, Michael Maire, Serge Belongie, James Hays, Pietro Perona, Deva Ramanan, Piotr Dollár, and C Lawrence Zitnick.  Microsoft coco: Common objects in context.  In _European conference on computer vision_ , pages 740–755. Springer, 2014.
  * [Ljosa et al.(2012)Ljosa, Sokolnicki, and Carpenter] Vebjorn Ljosa, Katherine L Sokolnicki, and Anne E Carpenter.  Annotated high-throughput microscopy image sets for validation.  _Nature methods_ , 9(7):637–637, 2012.
  * [Lu et al.(2015)Lu, Carneiro, and Bradley] Zhi Lu, Gustavo Carneiro, and Andrew P Bradley.  An improved joint optimization of multiple level set functions for the segmentation of overlapping cervical cells.  _IEEE Transactions on Image Processing_ , 24(4):1261–1272, 2015.
  * [Lu et al.(2016)Lu, Carneiro, Bradley, Ushizima, Nosrati, Bianchi, Carneiro, and Hamarneh] Zhi Lu, Gustavo Carneiro, Andrew P Bradley, Daniela Ushizima, Masoud S Nosrati, Andrea GC Bianchi, Claudia M Carneiro, and Ghassan Hamarneh.  Evaluation of three algorithms for the segmentation of overlapping cervical cells.  _IEEE journal of biomedical and health informatics_ , 21(2):441–450, 2016.
  * [Milletari et al.(2016)Milletari, Navab, and Ahmadi] Fausto Milletari, Nassir Navab, and Seyed-Ahmad Ahmadi.  V-net: Fully convolutional neural networks for volumetric medical image segmentation.  In _2016 fourth international conference on 3D vision (3DV)_ , pages 565–571. IEEE, 2016.
  * [Ren et al.(2015)Ren, He, Girshick, and Sun] Shaoqing Ren, Kaiming He, Ross Girshick, and Jian Sun.  Faster r-cnn: Towards real-time object detection with region proposal networks.  _Advances in neural information processing systems_ , 28, 2015.
  * [Ronneberger et al.(2015)Ronneberger, Fischer, and Brox] Olaf Ronneberger, Philipp Fischer, and Thomas Brox.  U-net: Convolutional networks for biomedical image segmentation.  In _International Conference on Medical image computing and computer-assisted intervention_ , pages 234–241. Springer, 2015.
  * [Schmidt et al.(2018)Schmidt, Weigert, Broaddus, and Myers] Uwe Schmidt, Martin Weigert, Coleman Broaddus, and Gene Myers.  Cell detection with star-convex polygons.  In _International Conference on Medical Image Computing and Computer-Assisted Intervention_ , pages 265–273. Springer, 2018.
  * [Wu et al.(2020)Wu, Chen, and Merhof] Yuli Wu, Long Chen, and Dorit Merhof.  Improving pixel embedding learning through intermediate distance regression supervision for instance segmentation.  In _ECCV Workshops_ , pages 213–227, 2020.
  * [Yu et al.()Yu, Lee, Hariharan, Bu, and Ahmed] Weimiao Yu, Hwee Kuan Lee, Srivats Hariharan, Wen Yu Bu, and Sohail Ahmed.  Ccdb: 6843, mus musculus, neuroblastoma.  _Cell Image Library._ URL https://doi.org/doi:10.7295/W9CCDB6843.