# Probabilistic degenerate Fubini polynomials associated with random variables

Rongrong Xu1,Taekyun Kim1,2,∗, Dae San Kim3,∗,and Yuankui Ma1,∗ 1School of
Science, Xi’an Technological University, Xi’an 710021, Shaanxi, China
2Department of Mathematics, Kwangwoon University, Seoul 01897, Republic of
Korea
3Department of Mathematics, Sogang University, Seoul 121-742, Republic of
Korea<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract.

Let $Y$ be a random variable such that the moment generating function of $Y$
exists in a neighborhood of the origin. The aim of this paper is to study
probabilistic versions of the degenerate Fubini polynomials and the degenerate
Fubini polynomials of order $r$, namely the probabilisitc degenerate Fubini
polynomials associated with $Y$ and the probabilistic degenerate Fubini
polynomials of order $r$ associated with $Y$. We derive some properties,
explicit expressions, certain identities and recurrence relations for those
polynomials. As special cases of $Y$, we treat the gamma random variable with
parameters $\alpha,\beta>0$, the Poisson random variable with parameter
$\alpha>0$, and the Bernoulli random variable with probability of success $p$.

###### Key words and phrases:

probabilistic degenerate Fubini polynomials; probabilistic degenerate Fubini
polynomials of order $r$

###### 2010 Mathematics Subject Classification:

11B73; 11B83

∗ Corresponding authors

##  1\. Introduction

In recent years, degenerate versions, $\lambda$-analogues and probabilistic
versions of many special polynomials and numbers have been investigated by
employing various methods such as generating functions, combinatorial methods,
umbral calculus, $p$-adic analysis, differential equations, probability,
special functions, analytic number theory and operator theory (see [11-16,
18-21] and the references therein).

Let $Y$ be a random variable satisfying the moment condition (see (17)). The
aim of this paper is to study probabilistic versions of the degenerate Fubini
polynomials and the degenerate Fubini polynomials of order $r$, namely the
probabilisitc degenerate Fubini polynomials associated with $Y$ and the
probabilistic degenerate Fubini polynomials of order $r$ associated with $Y$.
We derive some properties, explicit expressions, certain identities and
recurrence relations for those polynomials and numbers. In addition, we
consider the special cases that $Y$ is the gamma random variable with
parameters $\alpha,\beta>0$, the Poisson random variable with parameter
$\alpha(>0)$, and the Bernoulli random variable with probability of success
$p$.

The outline of this paper is as follows. In Section 1, we recall the
degenerate exponentials, the degenerate Stirling numbers of the second kind
${n\brace k}_{\lambda}$, the degenerate Bell polynomials, the degenerate
Fubini polynomials and the degenerate Fubini polynomials of order $r$. We
remind the reader of Lah numbers and the partial Bell polynomials. Assume that
$Y$ is a random variable such that the moment generating function of $Y$,
$E[e^{tY}]=\sum_{n=0}^{\infty}\frac{t^{n}}{n!}E[Y^{n}],\quad(|t|<r)$, exists
for some $r>0$. Let $(Y_{j})_{j\geq 1}$ be a sequence of mutually independent
copies of the random variable $Y$, and let
$S_{k}=Y_{1}+Y_{2}+\cdots+Y_{k},\,\,(k\geq 1)$, with $S_{0}=0$. Then we recall
the probabilistic degenerate Stirling numbers of the second kind associated
with $Y$ and the probabilistic degenerate Bell polynomials associated with
$Y$, $\phi_{n,\lambda}^{Y}(x)$. Also, we remind the reader of the gamma random
variable with parameters $\alpha,\beta>0$. Section 2 is the main result of
this paper. Let $(Y_{j})_{j\geq 1},\,\,S_{k},\,\,(k=0,1,\dots)$ be as in the
above. Then we first define the probabilistic degenerate Fubini polynomials
associated with the random variable $Y$, $F_{n,\lambda}^{Y}(x)$. We derive for
$F_{n,\lambda}^{Y}(x)$ an explicit expression in Theorem 2.1 and an expression
as an infinite sum involving $E[(S_{k})_{n,\lambda}]$ in Theorem 2.2. In
Theorem 2.3, when $Y\sim\Gamma(1,1)$, we find an expression for
$F_{n,\lambda}^{Y}(x)$ in terms of Lah numbers and Stirling numbers of the
first kind. We obtain a representation of $F_{n,\lambda}^{Y}(x)$ as an
integral over $(0,\infty)$ of the integrand involving
$\phi_{n,\lambda}^{Y}(x)$ and its generalization in Theorem 2.14. In Theorem
2.5, we express the probabilistic degenerate Fubini numbers associated with
$Y$, $F_{n,\lambda}^{Y}=F_{n,\lambda}^{Y}(1)$, as a finite sum involving the
partial Bell polynomials. Then we introduce the probabilistic degenerate
Fubini polynomials of order $r$ associated with $Y$ and deduce an explicit
expression for them in Theorem 2.6. We obtain a recurrence relation for
$F_{n,\lambda}^{Y}(x)$ in Theorem 2.7, and another one in Theorem 2.8 together
with its generalization in Theorem 2.15. In Theorem 2.9, the $r$th derivative
of $F_{n,\lambda}^{Y}(x)$ is expressed in terms of
$F_{i,\lambda}^{(r+1,Y)}(x)$. We get the identity
$\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{i=0}^{\infty}E[(S_{i})_{n,\lambda}]x^{i}$
in Theorem 2.10 and its generalization in Theorem 2.13. In Theorem 2.11, when
$Y$ is the Poisson random variable with parameter $\alpha$, we express
$F_{n,\lambda}^{Y}(x)$ in terms of the Fubini polynomials $F_{i}(x)$ and
${n\brace i}_{\lambda}$. In Theorem 12, when $Y$ is the Poisson random
variable with parameter $\alpha$, we show
$\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{k=0}^{\infty}\phi_{n,\lambda}(k\alpha)x^{k}$.
Finally, we show that $F_{n,\lambda}^{Y}(x)=F_{n,\lambda}(xp)$ if $Y$ is the
Bernoulli random variable with probability of success $p$. For the rest of
this section, we recall the facts that are needed throughout this paper.

For any $\lambda\in\mathbb{R}$, the degenerate exponentials are defined by

(1)
$e_{\lambda}^{x}(t)=\sum_{k=0}^{\infty}\frac{(x)_{k,\lambda}}{k!}t^{k},\quad
e_{\lambda}(t)=e_{\lambda}^{1}(t),\quad(\mathrm{see}\ [6-21]),$

where

(2)
$(x)_{0,\lambda}=1,\quad(x)_{n,.\lambda}=x(x-\lambda)\cdots\big{(}x-(n-1)\lambda\big{)},\quad(n\geq
1).$

Note that

$\lim_{\lambda\rightarrow 0}e_{\lambda}^{x}(t)=e^{xt}.$

The Stirling numbers of the first kind are defined by

(3) $(x)_{n}=\sum_{k=0}^{n}S_{1}(n,k)x^{k},\quad(n\geq 0),\quad(\mathrm{see}\
[1-3,5,24]),$

where

$(x)_{0}=1,\quad(x)_{n}=x(x-1)\cdots(x-n+1),\quad(n\geq 1).$

The Lah numbers are defined by

(4) $\langle x\rangle_{n}=\sum_{k=0}^{n}L(n,k)(x)_{k},\quad(n\geq
0),\quad(\mathrm{see}\ [5-24]),$

where

$\langle x\rangle_{0}=1,\quad\langle
x\rangle_{n}=x(x+1)\cdots(x+n-1),\quad(n\geq 1).$

By (4), we easily get

(5) $L(n,k)=\frac{n!}{k!}\binom{n-1}{k-1},\quad(n\geq k\geq
0),\quad(\mathrm{see}\ [5-24]).$

In [13], the degenerate Stirling numbers of the second kind are defined by

(6) $(x)_{n,\lambda}=\sum_{k=0}^{n}{n\brace k}_{\lambda}(x)_{k},\quad(n\geq
0).$

It is well known that the degenerate Bell polynomials are defined by

(7)
$e^{x(e_{\lambda}(t)-1)}=\sum_{n=0}^{\infty}\phi_{n,\lambda}(x)\frac{t^{n}}{n!},\quad(\mathrm{see}\
[12-14]).$

Thus, by (6) and (7), we get

(8) $\phi_{n,\lambda}(x)=\sum_{k=0}^{n}{n\brace k}_{\lambda}x^{k},\quad(n\geq
0),\quad(\mathrm{see}\ [12,17,21]).$

The degenerate Fubini polynomials are defined by

(9) $F_{n,\lambda}(x)=\sum_{k=0}^{n}{n\brace k}_{\lambda}k!x^{k},\quad(n\geq
0),\quad(\mathrm{see}\ [17-19,26]).$

Thus, by (9), we get

(10)
$\frac{1}{1-x(e_{\lambda}(t)-1)}=\sum_{n=0}^{\infty}F_{n,\lambda}(x)\frac{t^{n}}{n!},\quad(\mathrm{see}[4,19,21,26]).$

From (10), we note that

(11)
$\frac{1}{1-x}F_{n,\lambda}\bigg{(}\frac{x}{1-x}\bigg{)}=\bigg{(}x\frac{d}{dx}\bigg{)}_{n,\lambda}\frac{1}{1-x}=\sum_{k=0}^{\infty}(k)_{n,\lambda}x^{k},\quad(\mathrm{see}\
[18]).$

For $r\in\mathbb{N}$, the degenerate Fubini polynomials of order $r$ are
defined by

(12)
$\bigg{(}\frac{1}{1-y(e_{\lambda}(t)-1)}\bigg{)}^{r}=\sum_{n=0}^{\infty}F_{n,\lambda}^{(r)}(y)\frac{t^{n}}{n!},\quad(\mathrm{see}\
[8,9,16]).$

Thus, by (12), we get

(13) $F_{n,\lambda}^{(r)}(y)=\sum_{k=0}^{n}\binom{k+r-1}{k}y^{k}{n\brace
k}_{\lambda}k!,\quad(\mathrm{see}\ [8,18,19,22,23]).$

From (12), we have

(14)
$\bigg{(}\frac{1}{1-x}\bigg{)}^{r+1}F_{n,\lambda}^{(r+1)}\bigg{(}\frac{x}{1-x}\bigg{)}=\bigg{(}x\frac{d}{dx}\bigg{)}_{n,\lambda}\bigg{(}\frac{1}{1-x}\bigg{)}^{r+1}=\sum_{k=0}^{\infty}\binom{k+r}{r}(k)_{n,\lambda}x^{k},$

where $n,r$ are nonnegative integers.

For any integer $k\geq 0$, the partial Bell polynomials are given by

(15)
$\frac{1}{k!}\bigg{(}\sum_{i=1}^{\infty}x_{i}\frac{t^{i}}{i!}\bigg{)}^{k}=\sum_{n=k}^{\infty}B_{n,k}(x_{1},x_{2},\dots,x_{n-k+1})\frac{t^{n}}{n!},\quad(\mathrm{see}\
[5])$

where

(16) $\displaystyle B_{n,k}(x_{1},x_{2},\dots,x_{n-k+1})$
$\displaystyle=\sum_{\begin{subarray}{c}l_{1}+l_{2}+\cdots+l_{n-k+1}=k\\\
l_{1}+2l_{2}+\cdots+(n-k+1)l_{n-k+1}=n\end{subarray}}\frac{n!}{l_{1}l_{2}!\cdots
l_{n-k+1}!}\bigg{(}\frac{x_{1}}{1!}\bigg{)}^{l_{1}}\bigg{(}\frac{x_{2}}{2!}\bigg{)}^{l_{2}}\cdots\bigg{(}\frac{x_{n-k+1}}{(n-k+1)!}\bigg{)}^{l_{n-k+1}}.$

Let $Y$ be a random variable such that the moment generating function of $Y$

(17)
$E[e^{tY}]=\sum_{n=0}^{\infty}E[Y^{n}]\frac{t^{n}}{n!},\quad(|t|<r)\quad\textrm{exists
for some $r>0$}.$

Assume that $(Y_{j})_{j\geq 1}$ is a sequence of mutually independent copies
of $Y$ and $S_{k}=Y_{1}+Y_{2}+\cdots+Y_{k},\ (k\geq 1)$ with $S_{0}=0$.

The probabilistic degenerate Stirling numbers of the second kind associated
with random variable $Y$ are defined by

(18) ${n\brace
k}_{Y,\lambda}=\frac{1}{k!}\sum_{j=0}^{k}\binom{k}{j}(-1)^{k-j}E\big{[}(S_{j})_{n,\lambda}\big{]},\quad(n\geq
k\geq 0),\quad(\mathrm{see}\ [15]).$

By binomial inversion, the equation (18) is equivalent to

(19) $E\big{[}(S_{k})_{n,\lambda}\big{]}=\sum_{j=0}^{k}\binom{k}{j}j!{n\brace
j}_{Y,\lambda},\quad(\mathrm{see}\ [15]).$

From (18), we note that

(20)
$\frac{1}{k!}\Big{(}E[e_{\lambda}^{Y}(t)]-1\Big{)}^{k}=\sum_{n=k}^{\infty}{n\brace
k}_{Y,\lambda}\frac{t^{n}}{n!},\quad(k\geq 0),\quad(\mathrm{see}\ [15]).$

In view of (8), the probabilistic degenerate Bell polynomials associated with
$Y$ are defined by

(21) $\phi_{n,\lambda}^{Y}(x)=\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}x^{k},\quad(n\geq 0),\quad(\mathrm{see}\ [15,20]).$

When $Y=1$, we have $\phi_{n,\lambda}^{Y}(x)=\phi_{n,\lambda}(x)$.

By (21), we get

(22)
$e^{x(E[e_{\lambda}^{Y}(t)]-1}=\sum_{n=0}^{\infty}\phi_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!},\quad(\mathrm{see}\
[15]).$

We recall that $Y$ is the gamma random variable with parameter
$\alpha,\beta>0$ if probability density function of $Y$ is given by

$f(x)=\left\\{\begin{array}[]{cc}\frac{\beta}{\Gamma(\alpha)}e^{-\beta
x}(\beta x)^{\alpha-1},&\textrm{if $x\geq 0$},\\\ 0,&\textrm{if
$x<0$,}\end{array}\right.,\quad(\mathrm{see}\ [3,25,27]),$

which is denoted by $Y\sim\Gamma(\alpha,\beta)$.

##  2\. Probabilistic degenerate Fubini polynomials associated with random
variables

Let $(Y_{k})_{k\geq 1}$ be a sequence of mutually independent copies of random
variable $Y$, and let

$S_{0}=0,\quad S_{k}=Y_{1}+Y_{2}+\cdots+Y_{k},\quad(k\in\mathbb{N}).$

Now, we consider the probabilistic degenerate Fubini polynomials associated
with random variable $Y$ which are given by

(23)
$\frac{1}{1-x\big{(}E[e_{\lambda}^{Y}(t)]-1\big{)}}=\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}.$

For $Y=1$, we have $F_{n,\lambda}^{Y}(x)=F_{n,\lambda}(x),\ (n\geq 0)$. When
$x=1$, $F_{n,\lambda}^{Y}=F_{n,\lambda}^{Y}(1)$ are called the probabilistic
degenerate Fubini numbers associated with random variable $Y$.

From (23), we note that

(24)
$\displaystyle\frac{1}{1-x\big{(}E[e_{\lambda}^{Y}(t)]-1\big{)}}=\sum_{k=0}^{\infty}x^{k}\Big{(}E[e_{\lambda}^{Y}(t)]-1\Big{)}^{k}$
$\displaystyle=\sum_{k=0}^{\infty}x^{k}k!\frac{1}{k!}\Big{(}E\big{[}e_{\lambda}^{Y}(t)\big{]}-1\Big{)}^{k}=\sum_{k=0}^{\infty}x^{k}k!\sum_{n=k}^{\infty}{n\brace
k}_{Y,\lambda}\frac{t^{n}}{n!}$
$\displaystyle=\sum_{n=0}^{\infty}\sum_{k=0}^{n}x^{k}k!{n\brace
k}_{Y,\lambda}\frac{t^{n}}{n!}.$

Therefore, by comparing the coefficients on both sides of (24), we obtain the
following theorem.

######  Theorem 1.

For $n\geq 0$, we have

$F_{n,\lambda}^{Y}(x)=\sum_{k=0}^{n}{n\brace k}_{Y,\lambda}k!x^{k}.$

By (24), we get

(25)
$\displaystyle\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}=\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}=\frac{1}{1+x-xE[e_{\lambda}^{Y}(t)]}$
$\displaystyle=\frac{1}{1+x}\frac{1}{1-\frac{x}{1+x}E[e_{\lambda}^{Y}(t)]}=\frac{1}{1+x}\sum_{k=0}^{\infty}\bigg{(}\frac{x}{1+x}\bigg{)}^{k}\Big{(}E[e_{\lambda}^{Y}(t)]\Big{)}^{k}$
$\displaystyle=\sum_{n=0}^{\infty}\frac{1}{1+x}\sum_{k=0}^{\infty}\bigg{(}\frac{x}{1+x}\bigg{)}^{k}E\big{[}(Y_{1}+Y_{2}+\cdots+Y_{k})_{n,\lambda}\big{]}\frac{t^{n}}{n!}$
$\displaystyle=\sum_{n=0}^{\infty}\frac{1}{1+x}\sum_{k=0}^{\infty}\bigg{(}\frac{x}{1+x}\bigg{)}^{k}E[(S_{k})_{n,\lambda}]\frac{t^{n}}{n!}.$

Therefore, by (25), we obtain the following theorem.

######  Theorem 2.

For $n\geq 0$, we have

$F_{n,\lambda}^{Y}(x)=\frac{1}{1+x}\sum_{k=0}^{\infty}\bigg{(}\frac{x}{1+x}\bigg{)}^{k}E[(S_{k})_{n,\lambda}].$

In particular, for $Y=1$, we have

$F_{n,\lambda}^{Y}(x)=\frac{1}{1+x}\sum_{k=0}^{\infty}\bigg{(}\frac{x}{1+x}\bigg{)}^{k}(k)_{n,\lambda}.$

Let $Y\sim\Gamma(1,1)$. Then we have

(26)
$\displaystyle\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}=\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}=\sum_{k=0}^{\infty}x^{k}\Big{(}E[e_{\lambda}^{Y}(t)]-1\Big{)}^{k}$
$\displaystyle=\sum_{k=0}^{\infty}x^{k}\bigg{(}\int_{0}^{\infty}e_{\lambda}^{y}(t)e^{-y}dy-1\bigg{)}^{k}=\sum_{k=0}^{\infty}x^{k}\bigg{(}\int_{0}^{\infty}e^{y(\frac{1}{\lambda}\log(1+\lambda
t)-1)}dy-1\bigg{)}^{k}$
$\displaystyle=\sum_{k=0}^{\infty}k!x^{k}\frac{1}{k!}\bigg{(}\frac{\frac{1}{\lambda}\log(1+\lambda
t)}{1-\frac{1}{\lambda}\log(1+\lambda t)}\bigg{)}^{k}$
$\displaystyle=\sum_{k=0}^{\infty}k!x^{k}\sum_{l=k}^{\infty}L(l,k)\frac{1}{l!}\bigg{(}\frac{1}{\lambda}\log(1+\lambda
t)\bigg{)}^{l}$
$\displaystyle=\sum_{l=0}^{\infty}\sum_{k=0}^{l}k!x^{k}L(l,k)\sum_{n=l}^{\infty}\lambda^{n-l}S_{1}(n,l)\frac{t^{n}}{n!}$
$\displaystyle=\sum_{n=0}^{\infty}\sum_{l=0}^{n}\sum_{k=0}^{l}k!x^{k}L(l,k)\lambda^{n-l}S_{1}(n,l)\frac{t^{n}}{n!},$

where $S_{1}(n,l)$ are the Stirling numbers of the first kind. Here we should
observe that, for all $t$ with $|t|$ small, we have

$\big{|}\frac{1}{\lambda}\log(1+\lambda t)\big{|}<1,$

since $|\frac{\log(1+x)}{x}|$ is bounded on $(0,\infty)$. Therefore, by
comparing the coefficients on both sides of (26), we obtain the following
theorem.

######  Theorem 3.

Let $Y\sim\Gamma(1,1)$. Then we have

$F_{n,\lambda}^{Y}(x)=\sum_{l=0}^{n}\sum_{k=0}^{l}k!\lambda^{n-l}L(l,k)S_{1}(n,l)x^{k},\quad(n\geq
0).$

Now, we observe from (21) and Theorem 2.1 that

(27) $\displaystyle\int_{0}^{\infty}\phi_{n,\lambda}^{Y}(xy)e^{-y}dy$
$\displaystyle=\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}x^{k}\int_{0}^{\infty}y^{k}e^{-y}dy$
$\displaystyle=\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}x^{k}\Gamma(k+1)=\sum_{k=0}^{n}{n\brace k}_{Y,\lambda}x^{k}k!$
$\displaystyle=F_{n,\lambda}^{Y}(x),\quad(n\geq 0).$

Thus, from (27), we obtain the following theorem.

######  Theorem 4.

For $n\geq 0$, we have

$\int_{0}^{\infty}\phi_{n,\lambda}^{Y}(xy)e^{-y}dy=F_{n,\lambda}^{Y}(x).$

From (20), we note that

(28) $\displaystyle\sum_{n=k}^{\infty}{n\brace k}_{Y,\lambda}\frac{t^{n}}{n!}$
$\displaystyle=\frac{1}{k!}\Big{(}E[e_{\lambda}^{Y}(t)]-1\Big{)}^{k}=\frac{1}{k!}\bigg{(}\sum_{i=1}^{\infty}E[(Y)_{i,\lambda}]\frac{t^{i}}{i!}\bigg{)}^{k}$
$\displaystyle=\sum_{n=k}^{\infty}B_{n,k}\Big{(}E[(Y)_{1,\lambda}],E[(Y)_{2,\lambda}],\cdots,E[(Y)_{n-k+1,\lambda}]\Big{)}\frac{t^{n}}{n!}.$

Thus, by (28), we get

(29) ${n\brace
k}_{Y,\lambda}=B_{n,k}\Big{(}E[(Y)_{1,\lambda}],E[(Y)_{2,\lambda}],\cdots,E[(Y)_{n-k+1,\lambda}]\Big{)},\quad(n\geq
k\geq 0).$

Hence

(30) $\phi_{n,\lambda}^{Y}(y)=\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}y^{k}=\sum_{k=0}^{n}B_{n,k}\Big{(}E[(Y)_{1,\lambda}],E[(Y)_{2,\lambda}],\cdots,E[(Y)_{n-k+1,\lambda}]\Big{)}y^{k}.$

By (27) and (30), we get

(31) $\displaystyle F_{n,\lambda}^{Y}$
$\displaystyle=\int_{0}^{\infty}\sum_{k=0}^{n}B_{n,k}\Big{(}E[(Y)_{1,\lambda}],E[(Y)_{2,\lambda}],\cdots,E[(Y)_{n-k+1,\lambda}]\Big{)}\int_{0}^{\infty}y^{k}e^{-y}dy$
$\displaystyle=\sum_{k=0}^{n}k!B_{n,k}\Big{(}E[(Y)_{1,\lambda}],E[(Y)_{2,\lambda}],\cdots,E[(Y)_{n-k+1,\lambda}]\Big{)},\quad(n\geq
0).$

Therefore, by (31), we obtain the following theorem.

######  Theorem 5.

For $n\geq 0$, we have

$F_{n,\lambda}^{Y}=\sum_{k=0}^{n}k!B_{n,k}\Big{(}E[(Y)_{1,\lambda}],E[(Y)_{2,\lambda}],\cdots,E[(Y)_{n-k+1,\lambda}]\Big{)}.$

For $r\in\mathbb{N}$, the probabilistic degenerate Fubini polynomials of order
$r$ associated with random variable $Y$ are defined by

(32)
$\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}^{r}=\sum_{n=0}^{\infty}F_{n,\lambda}^{(r,Y)}(x)\frac{t^{n}}{n!}.$

When $Y=1$, we have $F_{n,\lambda}^{(r,Y)}(x)=F_{n,\lambda}^{(r)}(x),\ (n\geq
0)$.

From (20) and (32), we note that

(33)
$\displaystyle\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}^{r}=\sum_{i=0}^{\infty}\binom{-r}{i}(-1)^{i}x^{i}\Big{(}E[e_{\lambda}^{Y}(t)]-1\Big{)}^{i}$
$\displaystyle=\sum_{i=0}^{\infty}\binom{r+i-1}{i}i!x^{i}\frac{1}{i!}\Big{(}E[e_{\lambda}^{Y}(t)]-1\Big{)}^{i}=\sum_{i=0}^{\infty}\binom{r+i-1}{i}i!x^{i}\sum_{n=i}^{\infty}{n\brace
i}_{Y,\lambda}\frac{t^{n}}{n!}$
$\displaystyle=\sum_{n=0}^{\infty}\sum_{i=0}^{n}\binom{r+i-1}{i}i!x^{i}{n\brace
i}_{Y,\lambda}\frac{t^{n}}{n!}.$

Therefore, by (32) and (33), we obtain the following theorem.

######  Theorem 6.

For $n\geq 0$, we have

$F_{n,\lambda}^{(r,Y)}(x)=\sum_{i=0}^{n}\binom{r+i-1}{i}i!{n\brace
i}_{Y,\lambda}x^{i}.$

By (23), we get

(34)
$\displaystyle\sum_{n=1}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}=\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}-1=\frac{x(E[e_{\lambda}^{Y}(t)]-1)}{1-x(E[e_{\lambda}^{Y}(t)]-1)}$
$\displaystyle=\frac{xE[e_{\lambda}^{Y}(t)]}{1-x(E[e_{\lambda}^{Y}(t)]-1)}-\frac{x}{1-x(E[e_{\lambda}^{Y}(t)]-1)}$
$\displaystyle=x\sum_{k=0}^{\infty}E[(Y)_{k,\lambda}]\frac{t^{k}}{k!}\sum_{l=0}^{\infty}F_{l,\lambda}^{Y}(x)\frac{t^{l}}{l!}-x\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}.$
$\displaystyle=x\sum_{k=1}^{\infty}E[(Y)_{k,\lambda}]\frac{t^{k}}{k!}\sum_{l=0}^{\infty}F_{l,\lambda}^{Y}(x)\frac{t^{l}}{l!}$
$\displaystyle=\sum_{n=1}^{\infty}x\sum_{k=1}^{n}\binom{n}{k}E\big{[}(Y)_{k,\lambda}\big{]}F_{n-k,\lambda}^{Y}(x)\frac{t^{n}}{n!}.$

Therefore, by comparing the coefficients on both sides of (34), we obtain the
following theorem.

######  Theorem 7.

For $n\geq 1$, we have

$F_{n,}^{Y}(x)=x\sum_{k=1}^{n}\binom{n}{k}E[(Y)_{k,\lambda}]F_{n-k,\lambda}^{Y}(x).$

From (23), we note that

(35)
$\displaystyle\sum_{n=0}^{\infty}F_{n+1,\lambda}^{Y}(x)\frac{t^{n}}{n!}=\frac{d}{dt}\sum_{n=0}^{\infty}F_{n,,\lambda}^{Y}(x)\frac{t^{n}}{n!}=\frac{d}{dt}\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}$
$\displaystyle=\frac{xE[Ye_{\lambda}^{Y-\lambda}(t)]}{(1-x(E[e_{\lambda}^{Y}(t)]-1))^{2}}=\frac{x}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\frac{E[Ye_{\lambda}^{Y-\lambda}(t)]}{1-x(E[e_{\lambda}^{Y}(t)]-1)}$
$\displaystyle=x\sum_{i=0}^{\infty}F_{i,\lambda}^{Y}(x)\frac{t^{i}}{i!}\sum_{j=0}^{\infty}F_{j,\lambda}^{Y}(x)\frac{t^{j}}{j!}\sum_{m=0}^{\infty}E[Y(Y-\lambda)_{m,\lambda}]\frac{t^{m}}{m!}$
$\displaystyle=x\sum_{k=0}^{\infty}\bigg{(}\sum_{i=0}^{k}\binom{k}{i}F_{i,\lambda}^{Y}(x)F_{k-i,\lambda}^{Y}(x)\bigg{)}\frac{t^{k}}{k!}\sum_{m=0}^{\infty}E[(Y)_{m+1,\lambda}]\frac{t^{m}}{m!}$
$\displaystyle=\sum_{n=0}^{\infty}x\sum_{k=0}^{n}\sum_{i=0}^{k}\binom{k}{i}\binom{n}{k}F_{i,\lambda}^{Y}(x)F_{k-i,\lambda}^{Y}(x)E[(Y)_{n-k+1,\lambda}]\frac{t^{n}}{n!}.$

Therefore, by comparing the coefficients on both sides of (26), we obtain the
following theorem.

######  Theorem 8.

For $n\geq 0$, we have

$F_{n+1,\lambda}^{Y}(x)=x\sum_{k=0}^{n}\sum_{i=0}^{k}\binom{k}{i}\binom{n}{k}F_{i,\lambda}^{Y}(x)F_{k-i,\lambda}^{Y}(x)E[(Y)_{n-k+1,\lambda}].$

Now, we observe that

(36)
$\displaystyle\sum_{n=0}^{\infty}\frac{d^{r}}{dx^{r}}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}=\frac{d^{r}}{dx^{r}}\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}=r!\frac{\big{(}E[e_{\lambda}^{Y}(t)]\big{)}^{r}}{\big{(}1-x(E[e_{\lambda}^{Y}(t)]-1)\big{)}^{r+1}}$
$\displaystyle=r!\sum_{i=0}^{\infty}F_{i,\lambda}^{(r+1,Y)}(x)\frac{t^{i}}{i!}\sum_{k=0}^{\infty}E\big{[}(Y_{1}+Y_{2}+\cdots+Y_{r})_{k,\lambda}\big{]}\frac{t^{k}}{k!}$
$\displaystyle=\sum_{n=0}^{\infty}\bigg{(}r!\sum_{i=0}^{n}F_{i,\lambda}^{(r+1,Y)}(x)E[(S_{r})_{n-i,\lambda}]\binom{n}{i}\bigg{)}\frac{t^{n}}{n!}.$

Therefore, by (36), we obtain the following theorem.

######  Theorem 9.

For $r,n\geq 0$, we have

$\frac{d^{r}}{dx^{r}}F_{n,\lambda}^{Y}(x)=r!\sum_{i=0}^{n}F_{i,\lambda}^{(r+1,Y)}(x)E[(S_{r})_{n-i,\lambda}]\binom{n}{i}.$

From (19), we note that

(37)
$\displaystyle\sum_{n=0}^{\infty}\bigg{(}\sum_{i=0}^{\infty}E[(S_{i})_{n,\lambda}]x^{i}\bigg{)}\frac{t^{n}}{n!}=\sum_{i=0}^{\infty}x^{i}\sum_{n=0}^{\infty}E[(S_{i})_{n,\lambda}]\frac{t^{n}}{n!}$
$\displaystyle=\sum_{i=0}^{\infty}x^{i}E[e_{\lambda}^{S_{i}}(t)]=\sum_{i=0}^{\infty}x^{i}\Big{(}E[e_{\lambda}^{Y}(t)]\Big{)}^{i}$
$\displaystyle=\frac{1}{1-xE[e_{\lambda}^{Y}(t)]}=\frac{1}{1-x}\frac{1}{1-\frac{x}{1-x}\big{(}E[e_{\lambda}^{Y}(t)]-1\big{)}}$
$\displaystyle=\frac{1}{1-x}\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}\frac{t^{n}}{n!}.$

Therefore, by comparing the coefficients on both sides of (37), we obtain the
following theorem.

######  Theorem 10.

For $n\geq 0$, we have

$\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{i=0}^{\infty}E[(S_{i})_{n,\lambda}]x^{i}.$

Taking $x=\frac{1}{2}$, we get

$\sum_{i=0}^{\infty}E[(S_{i})_{n,\lambda}]\bigg{(}\frac{1}{2}\bigg{)}^{i}=2F_{n,\lambda}^{Y},\quad(n\geq
0).$

Let $Y$ be the Poisson random variable with parameter $\alpha(>0)$. Then we
have

(38)
$E[e_{\lambda}^{Y}(t)]=\sum_{n=0}^{\infty}e_{\lambda}^{n}(t)\frac{\alpha^{n}}{n!}e^{-\alpha}=e^{\alpha(e_{\lambda}(t)-1)}.$

From (38), we have

(39) $\displaystyle\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}$
$\displaystyle=\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}=\frac{1}{1-x(e^{\alpha(e_{\lambda}(t)-1)}-1)}$
$\displaystyle=\sum_{i=0}^{\infty}F_{i}(x)\alpha^{i}\frac{1}{i!}\big{(}e_{\lambda}(t)-1\big{)}^{i}=\sum_{i=0}^{\infty}F_{i}(x)\alpha^{i}\sum_{n=i}^{\infty}{n\brace
i}_{\lambda}\frac{t^{n}}{n!}$
$\displaystyle=\sum_{n=0}^{\infty}\bigg{(}\sum_{i=0}^{n}F_{i}(x){n\brace
i}_{\lambda}\alpha^{i}\bigg{)}\frac{t^{n}}{n!},$

where $F_{i}(x)$ are the Fubini polynomials given by
$\frac{1}{1-x(e^{t}-1)}=\sum_{i=0}^{\infty}F_{i}(x)\frac{t^{i}}{i!}$.
Therefore, by (39), we obtain the following theorem.

######  Theorem 11.

Let $Y$ be the Poisson random variable with parameter $\alpha(>0)$. Then we
have

$F_{n,\lambda}^{Y}(x)=\sum_{i=0}^{n}F_{i}(x){n\brace
i}_{\lambda}\alpha^{i},\quad(n\geq 0).$

Let $Y$ be the Poisson random variable with parameter $\alpha>0$. Then we have

(40)
$\Big{(}E[e_{\lambda}^{Y}(t)]\Big{)}^{k}=e^{k\alpha(e_{\lambda}(t)-1)}=\sum_{n=0}^{\infty}\phi_{n,\lambda}(k\alpha)\frac{t^{n}}{n!}.$

and

(41)
$\Big{(}E[e_{\lambda}^{Y}(t)]\Big{)}^{k}=\sum_{n=0}^{\infty}E[(S_{k})_{n,\lambda}]\frac{t^{n}}{n!}.$

Thus, by (40) and (41), we get

(42) $E\big{[}(S_{k})_{n,\lambda}\big{]}=\phi_{n,\lambda}(k\alpha),\quad(n\geq
0).$

From Theorem 2.10 and (42), we have

(43)
$\sum_{k=0}^{\infty}\phi_{n,\lambda}(k\alpha)x^{k}=\sum_{k=0}^{\infty}E\Big{[}(S_{k})_{n,\lambda}\Big{]}x^{k}=\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}.$

Therefore, by (43), we obtain the following theorem.

######  Theorem 12.

Let $Y$ be the Poisson random variable with parameter $\alpha(>0)$. For $n\geq
0$, we have

$\sum_{k=0}^{\infty}\phi_{n,\lambda}(k\alpha)x^{k}=\sum_{k=0}^{\infty}E\Big{[}(S_{k})_{n,\lambda}\Big{]}x^{k}=\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}.$

By using Theorem 2.6 and (19), we note that

(44)
$\displaystyle\bigg{(}\frac{1}{1-x}\bigg{)}^{r+1}F_{n,\lambda}^{(r+1,Y)}\bigg{(}\frac{x}{1-x}\bigg{)}$
$\displaystyle=\bigg{(}\frac{1}{1-x}\bigg{)}^{r+1}\sum_{l=0}^{n}{n\brace
l}_{Y,\lambda}\binom{l+r}{l}l!\bigg{(}\frac{x}{1-x}\bigg{)}^{l}$
$\displaystyle=\sum_{l=0}^{n}{n\brace
l}_{Y,\lambda}l!\binom{l+r}{l}x^{l}\bigg{(}\frac{1}{1-x}\bigg{)}^{l+r+1}$
$\displaystyle=\sum_{l=0}^{n}\sum_{k=0}^{\infty}{n\brace
l}_{Y,\lambda}l!\binom{l+r}{l}\binom{k+l+r}{k}x^{k+l}$
$\displaystyle=\sum_{l=0}^{n}\bigg{(}\sum_{k=l}^{\infty}{n\brace
l}_{Y,\lambda}l!\binom{l+r}{l}\binom{k+r}{k-l}\bigg{)}x^{k}$
$\displaystyle=\sum_{k=0}^{n}\binom{k+r}{k}x^{k}\sum_{l=0}^{k}{n\brace
l}_{Y,\lambda}(k)_{l}+\sum_{k=n+1}^{\infty}\binom{k+r}{k}x^{k}\sum_{l=0}^{n}{n\brace
l}_{Y,\lambda}(k)_{l}$
$\displaystyle=\sum_{k=0}^{n}\binom{k+r}{k}x^{k}\sum_{l=0}^{k}{n\brace
l}_{Y,\lambda}(k)_{l}+\sum_{k=n+1}^{\infty}\binom{k+r}{k}x^{k}\sum_{l=0}^{k}{n\brace
l}_{Y,\lambda}(k)_{l}$
$\displaystyle=\sum_{k=0}^{\infty}\binom{k+r}{k}x^{k}\sum_{l=0}^{k}{n\brace
l}_{Y,\lambda}(k)_{l}=\sum_{k=0}^{\infty}\binom{k+r}{k}x^{k}E\Big{[}(S_{k})_{n,\lambda}\Big{]}.$

Therefore, by (44), we obtain the following theorem.

######  Theorem 13.

For $n\geq 0$, we have

$\bigg{(}\frac{1}{1-x}\bigg{)}^{r+1}F_{n,\lambda}^{(r+1,Y)}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{k=0}^{\infty}\binom{k+r}{k}x^{k}E\Big{[}(S_{k})_{n,\lambda}\Big{]}.$

When $Y=1$, we have

$\bigg{(}\frac{1}{1-x}\bigg{)}^{r+1}F_{n,\lambda}^{(r+1)}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{k=0}^{\infty}\binom{k+r}{k}x^{k}(k)_{n,\lambda}.$

Now, we observe from (21) and Theorem 2.6 that

(45)
$\displaystyle\int_{0}^{\infty}y^{r-1}\phi_{n,\lambda}^{Y}(xy)e^{-y}dy=\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}x^{k}\int_{0}^{\infty}y^{r+k-1}e^{-y}dy$
$\displaystyle=\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}x^{k}\Gamma(r+k)=\Gamma(r)\sum_{k=0}^{n}{n\brace
k}_{Y,\lambda}x^{k}\binom{r+k-1}{k}k!$
$\displaystyle=\Gamma(r)F_{n,\lambda}^{(r,Y)}(x),\quad(r\in\mathbb{N}).$

Therefore, by (45), we obtain the following theorem.

######  Theorem 14.

For $n\geq 0$ and $r\geq 1$, we have

$F_{n,\lambda}^{(r,Y)}(x)=\frac{1}{\Gamma(r)}\int_{0}^{\infty}y^{r-1}\phi_{n,\lambda}^{Y}(xy)e^{-y}dy.$

From (32) and Theorem 2.14, we have

(46)
$\displaystyle\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}^{r}=\frac{1}{\Gamma(r)}\int_{0}^{\infty}y^{r-1}\sum_{n=0}^{\infty}\phi_{n,\lambda}^{Y}(xy)\frac{t^{n}}{n!}e^{-y}dy$
$\displaystyle=\frac{1}{\Gamma(r)}\int_{0}^{\infty}y^{r-1}e^{xy(E[e_{\lambda}^{Y}(t)]-1)}e^{-y}dy$
$\displaystyle=\frac{1}{\Gamma(r)}\int_{0}^{\infty}y^{r-1}e^{y(xE[e_{\lambda}^{Y}(t)]-1-x)}dy.$

By (46), we get

$\frac{1}{\Gamma(r)}\int_{0}^{\infty}y^{r-1}e^{y(xE[e_{\lambda}^{Y}(t)]-1-x)}dy=\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}^{r},$

where $r$ is a positive integer.

We observe that

(47) $\displaystyle\sum_{n=0}^{\infty}F_{n+1}^{(r,Y)}(x)\frac{t^{n}}{n!}$
$\displaystyle=\frac{d}{dt}\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}$
$\displaystyle=\frac{d}{dt}\bigg{(}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}\bigg{)}^{r}=\frac{rxE[e_{\lambda}^{Y-\lambda}(t)Y]}{\big{(}1-x(E[e_{\lambda}^{Y}(t)]-1)\big{)}^{r+1}}$
$\displaystyle=\frac{rx}{\big{(}1-x(E[e_{\lambda}^{Y}(t)]-1)\big{)}^{r}}\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}E[e_{\lambda}^{Y-\lambda}(t)Y]$
$\displaystyle=rx\sum_{l=0}^{\infty}F_{l,\lambda}^{(r,Y)}(x)\frac{t^{l}}{l!}\sum_{i=0}^{\infty}F_{i,\lambda}^{Y}(x)\frac{t^{i}}{i!}\sum_{j=0}^{\infty}E[(Y)_{j+1,\lambda}]\frac{t^{j}}{j!}$
$\displaystyle=rx\sum_{l=0}^{\infty}F_{l,\lambda}^{(r,Y)}(x)\frac{t^{l}}{l!}\sum_{k=0}^{\infty}\bigg{(}\sum_{j=0}^{k}\binom{k}{j}F_{k-j,\lambda}^{Y}(x)E[(Y)_{j+1,\lambda}]\bigg{)}\frac{t^{k}}{k!}$
$\displaystyle=\sum_{n=0}^{\infty}rx\sum_{k=0}^{n}\sum_{j=0}^{k}\binom{n}{k}\binom{k}{j}F_{n-k,\lambda}^{(r,Y)}(x)F_{k-j,\lambda}^{Y}(x)E[(Y)_{j+1,\lambda}]\frac{t^{n}}{n!}.$

Therefore, by (47), we obtain the following theorem.

######  Theorem 15.

For $n\geq 0$, we have

$F_{n+1,\lambda}^{(r,Y)}(x)=rx\sum_{k=0}^{n}\sum_{j=0}^{k}\binom{n}{k}\binom{k}{j}F_{n-k,\lambda}^{(r,Y)}(x)F_{k-j,\lambda}^{Y}(x)E[(Y)_{j+1,\lambda}].$

Let $Y$ be the Bernoulli random variable with probability of success $p$. Then
we have

(48)
$E[e_{\lambda}^{Y}(t)]=\sum_{k=0}^{1}e_{\lambda}^{k}(t)p(k)=1-p+pe_{\lambda}(t)=1+p(e_{\lambda}(t)-1).$

By (48), we get

(49) $\displaystyle\sum_{n=0}^{\infty}F_{n,\lambda}^{Y}(x)\frac{t^{n}}{n!}$
$\displaystyle=\frac{1}{1-x(E[e_{\lambda}^{Y}(t)]-1)}=\frac{1}{1-xp(e_{\lambda}(t)-1)}$
$\displaystyle=\sum_{n=0}^{\infty}F_{n,\lambda}(xp)\frac{t^{n}}{n!}.$

Therefore, by comparing the coefficients on both sides of (49), we obtain the
following theorem.

######  Theorem 16.

Let $Y$ be the Bernoulli random variable with probability of success $p$. For
$n\geq 0$, we have

$F_{n,\lambda}^{Y}(x)=F_{n,\lambda}(xp).$

##  3\. Conclusion

In this paper, we studied by using generating functions the probabilistic
degenerate Fubini polynomials associated with $Y$ and the probabilistic
degenerate Fubini polynomials of order $r$ associated with $Y$, as
probabilistic versions of the degenerate Fubini polynomials and the degenerate
Fubini polynomials of order $r$, respectively. Here $Y$ is a random variable
such that the moment generating function of $Y$ exists in a neighborhood of
the origin. In more detail, we derived several explicit expressions of
$F_{n,\lambda}^{Y}(x)$ (see Theorems 2.1, 2.2, 2.4) and those of
$F_{n,\lambda}^{r,Y}(x)$ (see Theorems 2.6, 2.14). We obtained a recurrence
relations for $F_{n,\lambda}^{Y}(x)$ (see Theorem 2.7), and another one (see
Theorem 2.8) together with its generalization (see Theorem 2.15). We expressed
the $r$th derivative of $F_{n,\lambda}^{Y}(x)$ in terms of
$F_{i,\lambda}^{(r+1,Y)}(x)$ (see Theorem 2.9). We showed the identity
$\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{i=0}^{\infty}E[(S_{i})_{n,\lambda}]x^{i}$
(see Theorem 2.10) and its generalization (see Theorem 2.13). We deduced an
explicit expression for $F_{n,\lambda}^{Y}(x)$ when $Y\sim\Gamma(1,1)$ (see
Theorem 2.3) and also that when $Y$ is the Poisson random variable with
parameter $\alpha$ (see Theorem 2.11). We proved that
$\frac{1}{1-x}F_{n,\lambda}^{Y}\bigg{(}\frac{x}{1-x}\bigg{)}=\sum_{k=0}^{\infty}\phi_{n,\lambda}(k\alpha)x^{k}$
when $Y$ is the Poisson random variable with parameter $\alpha$ (see Theorem
2.12). We showed $F_{n,\lambda}^{Y}(x)=F_{n,\lambda}(xp)$ when $Y$ be the
Bernoulli random variable with probability of success $p$ (see Theorem 2.16).

As one of our future projects, we would like to continue to study degenerate
versions, $\lambda$-analogues and probabilistic versions of many special
polynomials and numbers and to find their applications to physics, science and
engineering as well as to mathematics.

## Declarations

Ethical Statement The submitted article is an original research paper and has
not been published anywhere else. There are no clinical trials, animal
research, or human trials included in this paper.

Funding Tis research was funded by the National Natural Science Foundation of
China (No.12271320) and Key Research and Development Program of Shaanxi (No.
2023 ZDLGY-02).

Data Availability There is no data used in this article.

Conflict of Interest The authors declare no competing interests.

## References

  * [1] Abbas, M.; Bouroubi, S. _On new identities for Bell’s polynomials,_ Discrete Math. 293 (2005), no. 1-3, 5-10.
  * [2] Abramowitz, M.; Stegun, I. A. _Handbook of mathematical functions with formulas, graphs, and mathematical tables. For sale by the Superintendent of Documents,_ National Bureau of Standards Applied Mathematics Series, No. 55. U. S. Government Printing Office, Washington, DC, 1964. xiv+1046 pp.
  * [3] Adams, C. R.; Morse, A. P. _Random sampling in the evaluation of a Lebesgue integral,_ Bull. Amer. Math. Soc. 45 (1939), no. 6, 442-447.
  * [4] Adell, J. A.; Bényi, B.; Nkonkobe, S. _On higher order generalized geometric polynomials with shifted parameters,_ Quaest. Math. 46 (2023), no. 3, 551-567.
  * [5] Comtet, L. _Advanced combinatorics. The art of finite and infinite expansions._ Revised and enlarged edition. D. Reidel Publishing Co., Dordrecht, 1974. xi+343 pp. ISBN: 90-277-0441-4.
  * [6] Jang, G.-W.; Kim, T. _Some identities of central Fubini polynomials arising from nonlinear differential equation,_ Proc. Jangjeon Math. Soc. 22 (2019), no. 2, 205-213.
  * [7] Kilar, N.; Simsek, Y. _Identities and relations for Fubini type numbers and polynomials via generating functions and $p$-adic integral approach,_ Publ. Inst. Math. (Beograd) (N.S.) 106 (120) (2019), 113-123.
  * [8] Kim, D. S.; Jang, G.-W.; Kwon, H.-I.; Kim, T. _Two variable higher-order degenerate Fubini polynomials,_ Proc. Jangjeon Math. Soc. 21 (2018), no. 1, 5-22.
  * [9] Kim, D. S.; Kim, T.; Kwon, H.-I.; Park, J.-W. _Two variable higher-order Fubini polynomials,_ J. Korean Math. Soc. 55 (2018), no. 4, 975-986.
  * [10] Kim, D. S.; Kwon, J.; Dolgy, D. V.; Kim, T. _On central Fubini polynomials associated with central factorial numbers of the second kind,_ Proc. Jangjeon Math. Soc. 21 (2018), no. 4, 589-598.
  * [11] Kim, T.; Kim, D. S. _Some results on degenerate Fubini and degenerate Bell polynomials,_ Appl. Anal. Discrete Math. 17 (2023), no. 2, 548-560.
  * [12] Kim, T.; Kim, D. S. _Combinatorial identities involving degenerate harmonic and hyperharmonic numbers,_ Adv. in Appl. Math. 148 (2023), Paper No. 102535.
  * [13] Kim, T.; Kim,D. S. _Degenerate zero-truncated Poisson random variables,_ Russ. J. Math. Phys., 28 (2021), no. 1, 66-72.
  * [14] Kim, T.; Kim, D. S. _Some identities on truncated polynomials associated with degenerate Bell polynomials,_ Russ. J. Math. Phys., 28 (2021), no. 3, 342-355.
  * [15] Kim, T.; Kim, D. S. _Probabilistic degenerate Bell polynomials associated with random variables,_ Russ. J. Math. Phys. 30 (2023), no. 4, 528-542.
  * [16] Kim, T.; Kim, D. S.; Jang, G.-W. _A note on degenerate Fubini polynomials,_ Proc. Jangjeon Math. Soc. 20 (2017), no. 4, 521-531.
  * [17] Kim, T.; Kim, D. S.; Jang, G.-W.; Kim, D. _Two variable higher-order central Fubini polynomials,_ J. Inequal. Appl. 2019, Paper No. 146.
  * [18] Kim, T.; Kim, D. S.; Kim, H. K.; Lee, H. _Some properties on degenerate Fubini polynomials,_ Appl. Math. Sci. Eng. 30 (2022), no. 1, 235-248.
  * [19] Kim, T.; Kim, D. S.; Kwon, J. _Some identities related to degenerate $r$-Bell and degenerate Fubini polynomials,_ Appl. Math. Sci. Eng. 31 (2023), no. 1, Paper No. 2205642.
  * [20] Kim, T.; Kim, D. S.; Kwon, J. _Probabilistic degenerate stirling polynomials of the second kind and their applications,_ Math. Comput. Model. Dyn. Syst. (in press). Doi: 10.1080/13873954.2023.229757
  * [21] Kim, T.; Kim, D. S.; Lee, H.; Kwon, J. _On degenerate generalized Fubini polynomials,_ AIMS Math. 7 (2022), no. 7, 12227-12240.
  * [22] Kurt, B. _A parametric unified Apostol-type Bernoulli, Euler, Genocchi, Fubini polynomials and numbers,_ Filomat 37 (2023), no. 19, 6307-6317.
  * [23] Muhyi, A.; Araci, S. _A note on $q$-Fubini-Appell polynomials and related properties,_ J. Funct. Spaces 2021, Art. ID 3805809.
  * [24] Roman, S. _The umbral calculus,_ Pure and Applied Mathematics 111, Academic Press, Inc. [Harcourt Brace Jovanovich, Publishers], New York, 1984.
  * [25] Ross, S. M. _Introduction to probability models,_ Twelfth edition of Academic Press, London, 2019.
  * [26] Soni, R.; Pathak, A. K.; Vellaisamy, P. _A probabilistic extension of the Fubini polynomials,_ arXiv:2312.15193 [math.PR].
  * [27] Withers, C. S. _A simple expression for the multivariate Hermite polynomials,_ Statist. Probab. Lett. 47 (2000), no. 2, 165-169.


# BOD: Blindly Optimal Data Discovery

Thomas Hoang  Denison UniversityUSA<EMAIL_ADDRESS>
###### Abstract.

Combining discovery and augmentation is important in the era of data usage
when it comes to predicting the outcome of tasks. However, having to ask the
user the utility function to discover the goal to achieve the optimal small
rightful dataset is not an optimal solution. The existing solutions do not
make good use of this combination, hence underutilizing the data. In this
paper, we introduce a novel goal-oriented framework, called BOD: Blindly
Optimal Data Discovery, that involves humans in the loop and comparing utility
scores every time querying in the process without knowing the utility
function. This establishes the promise of using BOD: Blindly Optimal Data
Discovery for modern data science solutions.

PVLDB Reference Format:
PVLDB, 14(1): XXX-XXX, 2024.
doi:XX.XX/XXX.XX ††This work is licensed under the Creative Commons BY-NC-ND
4.0 International License. Visit https://creativecommons.org/licenses/by-nc-
nd/4.0/ to view a copy of this license. For any use beyond those covered by
this license, obtain permission by emailing<EMAIL_ADDRESS>Copyright is held
by the owner/author(s). Publication rights licensed to the VLDB Endowment.
Proceedings of the VLDB Endowment, Vol. 14, No. 1 ISSN 2150-8097.
doi:XX.XX/XXX.XX

PVLDB Artifact Availability:
The source code, data, and/or other artifacts have been made available at
https://github.com/Thomas12012002/BOD.

##  1\. Introduction

When data scientists have multiple datasets, each contains various variables,
and scientists want to know what they could get from the data, determining the
outcome after joining multiple relational tables. Assume all the datasets
augmenting are valid candidates, the values are integers, no overlapping
values with different attributes and each has no missing number value;
scientists often want to find the goal from the datasets after augmenting in
which they want to know what data can do with or predicting what based on the
datasets, for example housing prices are often based on location (which
contains near urban center or a place where the criminal status is low) and
goverment policies (which has taxes incentive) or the home values ( including
the age and size of the house). To predict the housing price, scientists with
their expertise in the field, rank each variable in each dataset one at a time
until there is no variable in any dataset remaining and augment the datasets
which are location, government policies in the area, home values.

The process of using machine learning application and causal relationship
tasks to predict the outcome requires a utility function that users have to
know exactly beforehand. This would be abundant for users to walk them through
the process of finding utility functions based on data in datasets when users
do not know their exact information. We solve this problem by introducing BOD:
Blindly Optimal Data Discovery by eliminating the finding utility function
process and by asking a small number of questions based on variables of
datasets to output a small subset of rightful tuples, which is the goal of
predicting the future.

The work of algorithms: The query first asks the user to pick the most
important variable in each dataset. For example, if a user want to predict
housing price without knowing their utility function, the query asks user to
rank each variable in each of the datasets once at a time in multiple rounds.
Based on the user’s ranking, the query selects all the tuples that meet the
user’s preferences based on utility score, and the process continues working
with these selected tuples. The query then outputs all the desirable tuples
whenever user stops ranking the variables.

Figure 1. Augmenting data from tables and filtering out items based on human
feedback of their ranking for variables

##  2\. Related Studies

Machine Learning. To make use of data to get to the point of predicting or
getting to know what the data can do, the technique of using machine learning
in the work (Galhotra et al., 2023) requires the user’s utility function.
However, this approach is not applicable in the case when the user’s utility
function is unknown. When in the case of an unknown utility function, the work
of (Gong et al., 2021) only selects the nearly optimal solution without
understanding what the data can do, our algorithm, the Blindly Optimal Data
Discovery, combines the advantages of these two works, which are getting to
know what the data can do and unknown user’s utility function.

Top-k. Top-K algorithm has many developments and ideas for a long time
(Fraternali et al., 2012) (Lee et al., 2009) (Lian and Chen, 2009) (Qin et
al., 2012) (Soliman et al., 2007), however, the algorithm requires user’s
utility function, which is not applicable in many cases and hard to obtain
itself.

Skyline. There are many studies related to Skyline algorithm including (Lee et
al., 2009) (Lin et al., 2007) (Mindolin and Chomicki, 2009) (Chan et al.,
2006a) (Chan et al., 2006b) (Tao et al., 2007) (Tao et al., 2009) (Papadias et
al., 2005) (Goncalves and Yidal, 2005) (Xia et al., 2008). The Skyline
algorithm, the disadvantage would be the output size, which we have no control
over, which means there may be uninterested tuples in the output and the
output may leave out the tuples that are dominated by others even the user
have interest in those. To fix this problem, researchers presented (Mouratidis
et al., 2021) combining the technique of top-k and Skyline together to control
the size of the output

Blindly Optimal Data Discovery. Our algorithm involves human interaction with
the data to rank the variables in each of the datasets once at a time, and we
filter out the promising tuples considering the utility score after each
ranking phase.

##  3\. Problem Definition

Table 1. Notation and meaning $D$ | The set of all tuples
---|---
$d$ | The total attributes of the final table after augmenting all the datasets
$S$ | The table of all the relational tables augmented over $n_{i}$ attributes
$t$ | The specific tuple that has an instance schema
$y_{min_{i}}$ | The maximum values’s summation based on user’s $i$th ranking
$y_{max_{i}}$ | The summation of all values based on user’s $i$th ranking
$y_{i}$ | The summation of all values of the instance

Let $(S_{1},\dots,S_{n})$ denote a relation table over $n$ attributes, $S_{i}$
denotes the $i$th attribute. $D$ composes the list of tuples with a table
$(S_{1},\dots,S_{n})$ such that every tuples $t\in T$ is an instance of the
table.

Linear function:

  * •

$LINEAR=\\{f|f(x)=\sum\limits_{i=1}^{d}f(x_{i})\\\ \text{\\{where each}\
f(x_{i})\text{ is a linear function.}\\}$

Assume we want to predict the housing price in an area, we have datasets
including location, government policies in the area, and home values. In
detail, the dataset location composes a relational table that has two
attributes that are Near Urban and Criminal Free, the dataset government
policies composes table that has one attribute which is Tax, while the dataset
home values composes a tables that has two attributes which are Size and Age
of the house. After augmenting all the datasets, we have a dataset composes
five attributes and six instances after augmenting the tables as shown in the
table 2

Price | Near Urban | Criminal Free | Tax | Size | Age
---|---|---|---|---|---
House 1 | 26 | 5 | 90 | 1500 | 150
house 2 | 35 | 2 | 20 | 1300 | 120
House 3 | 45 | 3 | 78 | 2000 | 95
House 4 | 9 | 2 | 46 | 1700 | 50
House 5 | 6 | 4 | 65 | 1800 | 25
house 6 | 47 | 7 | 30 | 1450 | 75
Table 2. Example of the tables after augmenting

Then we scale down the values in the range of $\\{0,1\\}$ by the highest value
in each variables and we calculate $y_{min_{i}}$ after asking the user to rank
the variables in the original tables as shown in table 3.

Price | Near Urban | Criminal Free | Tax | Size | Age | Sum
---|---|---|---|---|---|---
House 1 | 0.55 | 0.71 | 1 | 0.75 | 1 | 4.01
house 2 | 0.74 | 0.29 | 0.22 | 0.65 | 0.8 | 2.7
House 3 | 0.96 | 0.43 | 0.87 | 1 | 0.63 | 3.89
House 4 | 0.19 | 0.29 | 0.51 | 0.85 | 0.33 | 2.17
House 5 | 0.13 | 0.57 | 0.72 | 0.9 | 0.17 | 2.49
house 6 | 1 | 1 | 0.33 | 0.73 | 0.5 | 3.56
Table 3. Example of the query scaling down all the values in all variables
after augmenting

Suppose the user rank Near urban , Tax, and Size as the first priority in the
first interaction with the algorithm, we then calculate the summation of each
instance and find the maximum calculation in this round, which is
$y_{min_{1}}=2.83$ in this case as shown in table 4.

Price | Near Urban | Tax | Size | Near Urban + Tax + Size
---|---|---|---|---
House 1 | 0.55 | 1 | 0.75 | 2.3
house 2 | 0.74 | 0.22 | 0.65 | 1.61
House 3 | 0.96 | 0.87 | 1 | 2.83
House 4 | 0.19 | 0.51 | 0.85 | 1.55
House 5 | 0.13 | 0.72 | 0.9 | 1.75
house 6 | 1 | 0.33 | 0.73 | 2.06
Table 4. Example of the query when user rank Near Urban, Tax, and Size as the
first priorities

Based on the information from the first interaction, we find the tuples that
have $y_{min_{1}}\leq y_{1}\leq y_{max_{1}}$ that are House 3 and House 6
since $2.83\leq y_{1}\leq 3.89$. The table 5 shown the result from the first
interaction. The final interaction which is the second ranking would result in
only House 3.

Price | Near Urban | Criminal Free | Tax | Size | Age | Sum
---|---|---|---|---|---|---
House 3 | 0.96 | 0.43 | 0.87 | 1 | 0.63 | 3.89
house 6 | 1 | 1 | 0.33 | 0.73 | 0.5 | 3.56
Table 5. The results of the algorithm

##  4\. BOD: Blindly Optimal Data Discovery

Since we scale down each value in a variable by a factor of max value of the
variable such that every value now in the variable is in the range of
$\\{0,1\\}$ and we scale the values in every variable after augmenting all the
datasets. For each of the relational tables, We ask the user to pick the first
most important variable in each of dataset, knowing that the number of
variables in each table is varied, called $n_{i}$. We then find the maximum
sum of values from the most first important variables from the original tables
after augmenting, called $y_{min_{1}}$. Then, we find $y_{max_{1}}$ from the
row with the maximum sum $y_{min_{1}}$, which is the sum of all values in that
row. After that, we select the items that have their utility score, $y_{1}$ ,
which is the sum of all variables, in which $y_{min_{1}}\leq y_{1}\leq
y_{max_{1}}$. We then use these items in a new dataset and begin the second
phase, and so on, until we get to the last variable in any table or the user
wants to stop ranking the variables in each dataset.

Algorithm 1 Blindly Optimal Data Discovery

0: The tables $S_{i}$ that has $D$ tuples.

0: return a subset $t$ of $D$.

1: Let $S_{i}=1$

2: while $S_{i}\leq$ length of $S$ do

3: for $i$ in the length of $S_{i}$ do

4: Let the user rank the attributes in each table $S_{i}$ until the user stops

5: $S_{i}+=1$

6: Find the maximum length of each $S_{i}$ in $S$, called $max_{n_{i}}$

7: Concatenate DataFrames/Tables horizontally along columns

8: $d_{i}=1$

9: for $d_{i}\leq d$ do

10: find max value called $S_{m}ax_{i}$

11: Scale down all values by dividing $S_{max_{i}}$

12: $d_{i}+=1$

13: Update the new $S$

14: for $i$ in range $max_{n_{i}}$ do

15: Based one user ranking in $i$th time, calculate $y_{min_{i}}$ and
$y_{max_{i}}$ accordingly

16: If $y_{min_{i}}\leq y_{i}\leq y_{max_{i}}$

17: Select all tuples $t$ that satisfies and let update $S$

18: return $t$

##  5\. BOD’s Performance Guarantee

We show that Blindly Optimal Data Discovery can identify the optimal solution
in exactly $O(d)$ queries and returns a nearly optimal solutions in a maximum
of $O(d)$ queries since the number of questions depends on the amount of the
total attributes after augmenting all the datasets. The performance runs far
better than any previously known algorithms.

##  6\. Experimental Evaluation

We assume the augmentation part is correct with no error from the scientist’s
domain; all values are valid, not missing, and all are in positive integers.
We test the algorithm using Google Collaboration running on GL65 Leopard
10SCXK, an x64-based PC, on Microsoft Windows 11 Home Single Language. We test
on three synthetic datasets, each dataset has 10000 tuples and the range of 3
to 9 variables. The values in each cell are randomly generated in the range of
$\\{1,1000\\}$.

Figure 2. Changes in number of attributes

Then, We test on three synthetic datasets in the range of 5000, 10000, and
15000 tuples, each with six variables. The values in each cell are randomly
generated in the range of $\\{1,1000\\}$.

Figure 3. Changes in number of tuples

##  7\. Conclusion

The advantage is that the algorithm does not depend on the utility function
and asks users the maximum amount of questions based on the total attributes
after augmenting all the datasets, which is $d$. The algorithm outputs a small
subset of items that help predict the optimal outcomes scientists could use
with data science applications. The running time that we test in the case of
three synthetic datasets, each dataset has 10000 tuples and a range of 3 to 9
variables. The values in each cell are randomly generated in the range of
$\\{1,1000\\}$ shown in Figure 2, have the maximum running time under 30
seconds while the same maximum running time happens in the case of three
synthetic datasets in the range of 5000, 10000, and 15000 tuples, each with
six variables. The values in each cell are randomly generated in the range of
$\\{1,1000\\}$.

##  8\. Citations

## References

  * (1)
  * Chan et al. (2006a) C. Chan, H. Jagadish, K. Tan, A. Tung, and Z. Zhang. 2006a.  Finding k-dominant skylines in high dimensional space. In _Proceedings of the ACM SIGMOD International Conference on Management of Data_.
  * Chan et al. (2006b) C. Chan, H. Jagadish, K. Tan, A. Tung, and Z. Zhang. 2006b.  On high dimensional skylines. In _Advances in Database Technology-EDBT 2006_. Springer, 478–495.
  * Fraternali et al. (2012) P. Fraternali, D. Martinenghi, and M. Tagliasacchi. 2012.  Top-k bounded diversification. In _Proceedings of the 2012 International Conference on Management of Data_.
  * Galhotra et al. (2023) Sainyam Galhotra, Yue Gong, and Raul Castro Fernandez. 2023.  Metam: Goal-Oriented Data Discovery.  _2023 IEEE 39th International Conference on Data Engineering (ICDE)_ (2023), 2780–2793.  https://api.semanticscholar.org/CorpusID:258187398
  * Goncalves and Yidal (2005) M. Goncalves and M. Yidal. 2005.  Top-k skyline: a unified approach. In _On the Move to Meaningful Internet System 2005_.
  * Gong et al. (2021) Yue Gong, Zhiru Zhu, Sainyam Galhotra, and Raul Castro Fernandez. 2021.  Ver: View Discovery in the Wild.  _2023 IEEE 39th International Conference on Data Engineering (ICDE)_ (2021), 503–516.  https://api.semanticscholar.org/CorpusID:252545411
  * Lee et al. (2009) Jongwuk Lee, Gae won You, and Seung won Hwang. 2009.  Personalized top-k skyline queries in high-dimensional space.  _Information Systems_ 34, 1 (2009), 45–61.
  * Lian and Chen (2009) X. Lian and L. Chen. 2009.  Top-k dominating queries in uncertain databases. In _Proceedings of International Conference on Extending Database Technology: Advances in Database Technology_.
  * Lin et al. (2007) X. Lin, Y. Yuan, Q. Zhang, and Y. Zhang. 2007.  Selecting stars: The k most representative skyline operator. In _Proceedings of International Conference on Data Engineering_.
  * Mindolin and Chomicki (2009) D. Mindolin and J. Chomicki. 2009.  Discovering relative importance of skyline attributes. In _Proceedings of the VLDB Endowment_.
  * Mouratidis et al. (2021) Kyriakos Mouratidis, Keming Li, and Bo Tang. 2021.  Marrying Top-k with Skyline Queries: Relaxing the Preference Input While Producing Output of Controllable Size. In _Proceedings of the 2021 International Conference on Management of Data_ (Virtual Event, China) _(SIGMOD/PODS ’21)_. Association for Computing Machinery, New York, NY, USA, 1317–1330.
  * Papadias et al. (2005) D. Papadias, Y. Tao, G. Fu, and B. Seeger. 2005.  Progressive skyline computation in database systems. In _ACM Transactions on Database Systems (TODS)_ , Vol. 30. ACM, 41–82.
  * Qin et al. (2012) L. Qin, J. Yu, and L. Chang. 2012.  Diversifying top-k results. In _Proceedings of the VLDB Endowment_.
  * Soliman et al. (2007) M. Soliman, I. Ilyas, and K. Chen-Chuan Chang. 2007.  Top-k query processing in uncertain databases. In _Proceedings of International Conference on Data Engineering_. IEEE, 896–905.
  * Tao et al. (2009) Y. Tao, L. Ding, and J. Pei. 2009.  Distance-based representative skyline. In _Proceedings of International Conference on Data Engineering_.
  * Tao et al. (2007) Y. Tao, X. Xiao, and J. Pei. 2007.  Efficient Skyline and Top-k Retrieval in Subspaces. In _TKDE_.
  * Xia et al. (2008) T. Xia, D. Zhang, and Y. Tao. 2008.  On skylineing with flexible dominance relation. In _Proceedings of International Conference on Data Engineering_.


# Spontaneous symmetry breaking-induced thermospin effect in superconducting
tunnel junctions

Gaia Germanese  Dipartimento di Fisica dell’Università di Pisa, Largo
Pontecorvo 3, I-56127 Pisa, Italy  NEST Istituto Nanoscienze-CNR and Scuola
Normale Superiore, I-56127 Pisa, Italy  Federico Paolucci  NEST Istituto
Nanoscienze-CNR and Scuola Normale Superiore, I-56127 Pisa, Italy  Giampiero
Marchegiani  Quantum Research Centre, Technology Innovation Institute, Abu
Dhabi, UAE  Alessandro Braggio  NEST Istituto Nanoscienze-CNR and Scuola
Normale Superiore, I-56127 Pisa, Italy  Francesco Giazotto  NEST Istituto
Nanoscienze-CNR and Scuola Normale Superiore, I-56127 Pisa, Italy

###### Abstract

We discuss the charge and the spin tunneling currents between two Bardeen-
Cooper-Schrieffer (BCS) superconductors, where one density of states is spin-
split. In the presence of a large temperature bias across the junction, we
predict the generation of a spin-polarized thermoelectric current. This
thermo-spin effect is the result of a spontaneous particle-hole symmetry
breaking in the absence of a polarizing tunnel barrier. The two spin
components, which move in opposite directions, generate a spin current larger
than the purely polarized case when the thermo-active component dominates over
the dissipative one.

Figure 1: (a) Top: $S-I-S_{m}$ junction building blocks. The left hot [right
cold] superconductor is indicated with $L$ ($R$), and the exchange field is
present in $S_{m}$ only. Bottom: energy band diagrams of the superconductors
are shown in the absence ($h_{exc}=0$, left panel) and presence
($h_{exc}=0.4\Delta_{0}$, right panel) of an exchange field ($h_{exc}$). For
$h_{exc}=0$, the DoS of R exhibits two peaks. For $h_{exc}=0.4\Delta_{0}$, the
spin-split DoS (blue $\uparrow$ and green $\downarrow$ for the up and down
spin components, respectively) shows four peaks. The difference between the
chemical potential of the two leads is $\mu_{R}-\mu_{L}=eV$. The temperature
of the $L$ ($R$) superconductor is set at $T_{L}=0.7T_{c,L}$
($T_{R}=0.01T_{c,L}$). (b) Quasiparticle current-voltage characteristics for
$h_{exc}=0$ (black) and $h_{exc}=0.4\Delta_{0}$ (aquamarine) display two and
four matching peaks, respectively. Inset: Blow-up of the charge ($I_{q}$,
solid lines) and spin currents ($I_{s}$, dashed lines) at low bias voltage.

Introduction. Spintronics consists of the active manipulation of the spin
degree of freedom to develop a wide range of solid-state technologies (Žutić
_et al._ , 2004; Linder and Robinson, 2015). In this context, hybrid systems
involving ferromagnets and superconductors have been exploited for the
generation of spin-polarized currents (Tedrow and Meservey, 1971; Meservey
_et al._ , 1975; Takahashi  _et al._ , 1999; Jedema  _et al._ , 1999;
Takahashi  _et al._ , 2001; Yamashita  _et al._ , 2002; Maekawa  _et al._ ,
2002; Giazotto and Taddei, 2008) with applications as memory elements (Hao
_et al._ , 1990; Giazotto, 2009; Bergeret  _et al._ , 2012; De Simoni  _et
al._ , 2018), due to a nearly perfect spin-valve effect. At the same time,
thermoelectricity in these hybrid devices has been discussed theoretically
(Virtanen and Heikkilä, 2004; Titov, 2008; Jacquod and Whitney, 2010; Kalenkov
and Zaikin, 2017; Heidrich and Beckmann, 2019) and experimentally (Eom  _et
al._ , 1998; Parsons  _et al._ , 2003; Jiang and Chandrasekhar, 2005; Giazotto
_et al._ , 2015; Tan  _et al._ , 2021) with intriguing results at non-local
level (Hussein  _et al._ , 2019; Kirsanov  _et al._ , 2019; Blasi  _et al._ ,
2020a, b) and applications as detectors Heikkilä  _et al._ (2018). Yet, in the
presence of temperature gradients, the interplay of the magnetic field with
the superconducting order parameter gives rise to exotic non-equilibrium
phenomena, such as the generation of pure spin currents (spin-Seebeck effect)
(Machon  _et al._ , 2013; Ozaeta  _et al._ , 2014; Linder and Bathen, 2016;
Bergeret  _et al._ , 2018; Kolenda  _et al._ , 2017; Shapiro  _et al._ ,
2017).
Particle-hole (PH) symmetry breaking is a necessary requirement to generate
thermoelectric power in the linear regime, i.e., for a small voltage ($\delta
V$) and temperature ($\delta T$) bias. This condition results in the following
inequality for the charge current, $I_{q}(\delta V,\delta T)\neq-I_{q}(-\delta
V,\delta T)$ Sánchez and López (2016); Benenti  _et al._ (2017). This
limitation can be overcome in the presence of a large temperature bias
(nonlinear regime), as recently demonstrated in tunnel junctions between
superconductors with different energy gaps (Marchegiani  _et al._ , 2020a, b).
Here, we investigate the thermo-spin effect induced by spontaneous PH symmetry
breaking in hybrid ferromagnetic insulator/superconductor systems for a large
temperature bias. In these conditions, dissipative and thermo-active opposite
spin components coexist, resulting in maximal spin current efficiency and
thermoelectric power generation. We consider a thermally-biased tunnel
junction between a Bardeen-Cooper-Schrieffer (BCS) superconductor ($S$) and a
ferromagnetic superconductor ($S_{m}$), labeled as left ($L$) and right ($R$),
respectively [see Fig. 1(a)]. The $S_{m}$ electrode can be realized by the
proximity of a ferromagnetic insulator, which induces an exchange interaction
(Tokuyasu  _et al._ , 1988; Hao  _et al._ , 1991; Strambini  _et al._ , 2017),
or, eventually, by an in-plane magnetic field that induces Zeeman splitting
without suppressing superconductivity in the film (Tedrow  _et al._ , 1986;
Moodera  _et al._ , 1988; Tinkham, 2004). In this system, the exchange
interaction ($h_{exc}$) resolves the degeneracy between the spin up
($\uparrow$) and spin down ($\downarrow$) components of the superconducting
Density of States (DoS) 111We assume the insulating barrier is sufficiently
opaque suppressing the Josephson coupling or the in-plane magnetic field
suppresses it by Franhoufer interference Marchegiani  _et al._ (2020c); Barone
(1992)..
We first give an intuitive description of the effect by analyzing the charge
current ($I_{q}$) and the spin current ($I_{s}$) in the presence of a thermal
bias, i.e., $T_{L}>T_{R}$, where $T_{L}~{}(T_{R})$ is the temperature of the
left (right) electrode. Figure 1(b) shows $I_{q}$ (solid curves) and $I_{s}$
(dashed curves) as a function of the bias voltage ($V$) in the absence [black,
see also left energy diagram in the box of Fig. 1(a)], and in the presence
[aquamarine, see also right energy diagram in the box of Fig. 1(a)] of a
sizeable $h_{\rm exc}$. In both cases, $I_{q}$ is linear in the voltage bias
($I_{q}\propto V$) for large $V$ (ohmic response). For low values of $V$, the
trace is highly nonlinear and current peaks appear, corresponding to the
matching between the singularities in the superconducting DoS Shapiro  _et
al._ (1962). Note that each of the two antisymmetric peaks in the $I_{q}(V)$
characteristics present for $h_{\rm exc}=0$ is doubled for $h_{\rm exc}\neq
0$, since the spin-degeneracy is broken by the exchange interaction. Indeed,
$I_{s}$ is exactly zero for $h_{\rm exc}=0$, and finite for $h_{\rm exc}\neq
0$. Moreover, $I_{q}$ flows against the bias voltage for small values of $V$
when $h_{\rm exc}\neq 0$, i.e., the junction produces thermoelectric power,
while it is always dissipative for $h_{\rm exc}=0$. When the active thermo-
spin current is generated, it results $|I_{s}|>|I_{q}|$. This inequality is
never realized in a perfect spin-polarized barrier (i.e. $100\%$
polarization), where $|I_{s}|=|I_{q}|$ Žutić  _et al._ (2004). Indeed, in our
structure, the two finite spin-current components ($\rm{I_{\uparrow}}$ and
$\rm{I_{\downarrow}}$) flow in opposite directions. Moreover, when the thermo-
active component dominates, a net thermoelectric spin-polarized current
occurs.
Model. We assume $S_{m}$ to be much thinner than its superconducting coherence
length ($\xi_{0,R}$), in order to consider a homogeneous exchange interaction
in $S_{m}$ Tokuyasu  _et al._ (1988). The normalized DoS of the spin-$\sigma$
component [with $\sigma=+(\uparrow),-(\downarrow)$] for the $\alpha$
superconducting electrode (with $\alpha=L,R$) is expressed as

$N_{\alpha\sigma}(E_{\alpha})=\left|\Re\\!\\!\left[\frac{E_{\alpha}+\sigma
h_{\alpha}+i\Gamma}{\sqrt{(E_{\alpha}+\sigma
h_{\alpha}+i\Gamma)^{2}-{\Delta_{\alpha}}^{2}}}\right]\\!\\!\right|,$ (1)

where $h_{\alpha}$ is the exchange field, $E_{\alpha}=E-\mu_{\alpha}$ is the
quasiparticle energy measured with respect to the chemical potential
$\mu_{\alpha}$ 222For simplicity, we assume that the chemical potential for
the two spin-species is equal $\mu_{\uparrow}=\mu_{\downarrow}$ (no spin
accumulation). This approximation is valid for sufficiently opaque tunneling
barrier., $\Delta_{\alpha}$ is the self-consistent superconducting energy gap
for the $\alpha$-lead 333The self-consistency relation is Giazotto (2009);
Giazotto  _et al._ (2015):
$\log\left(\frac{\Delta_{0}}{\Delta_{\alpha}}\right)=\int_{0}^{\hbar\omega_{D}}dE\,\frac{f_{+}(E)-f_{-}(E)}{\sqrt{E^{2}+\Delta_{\alpha}^{2}}},$
(2) where $\Delta_{0}$ is the zero-temperature pairing potential in absence of
the exchange field, $\omega_{D}$ is the Debye frequency and
$f_{\pm}(E)=(1+\exp[(\sqrt{E^{2}+\Delta^{2}}\mp h_{\rm
exc})/k_{B}T_{R}])^{-1}$ with $k_{B}$ Boltzmann constant., and $\Gamma$ is the
phenomenological Dynes parameter Dynes  _et al._ (1984). The sum of two spin
contributions gives the total DoS of a spin-split superconductor
$N_{\alpha}(E)=\sum_{\sigma=\pm}N_{\alpha\sigma}(E)/2$. For $\alpha=L$, we
assume $h_{L}=0$. Note that, even in the presence of an exchange field, the
DoS is PH symmetric, i.e., $N_{\alpha}(E_{\alpha})=N_{\alpha}(-E_{\alpha})$
Tinkham (2004), and spin component satisfies
$N_{\alpha\sigma}(E_{\alpha})=N_{\alpha\bar{\sigma}}(-E_{\alpha})$, with
$\bar{\sigma}=-\sigma$.
The quasiparticle current of the spin-$\sigma$ component reads

$I_{\sigma}=\frac{G_{T}}{2e}\\!\\!\int_{-\infty}^{\infty}\\!\\!\\!\\!\\!\\!\\!\\!dE~{}N_{L\sigma}(E-eV)N_{R\sigma}(E)F_{LR}(E),$
(3)

where $e$ is the electron charge, $G_{T}$ is the normal-state conductance of
the tunnel junction 444$G_{T}$ is defined as
$G_{T}=G_{\uparrow}+G_{\downarrow}$. We assume that the two spin contributes
are equal, such as $G_{\uparrow}=G_{\downarrow}$., and
$F_{LR}(E)=f_{L}(E-eV)-f_{R}(E)$ is the difference between the Fermi-Dirac
quasiparticle distributions of two electrodes
$f_{\alpha}(E,T_{\alpha})=[1+\exp(E/k_{B}T_{\alpha})]^{-1}$ with $\alpha=L,R$.
The charge current and the spin current are defined as
$I_{q}=I_{\uparrow}+I_{\downarrow}$ and $I_{s}=I_{\uparrow}-I_{\downarrow}$,
respectively. By exploiting the PH symmetry, we note that the charge current
is odd (even) in $V$ ($h_{\rm exc}$) with $I_{q}(V,h_{\rm
exc})=-I_{q}(-V,h_{\rm exc})=I_{q}(V,-h_{\rm exc})$, while the spin current is
an even (odd) function in $V$ ($h_{\rm exc}$) with $I_{s}(V,h_{\rm
exc})=I_{s}(-V,h_{\rm exc})=-I_{s}(V,-h_{\rm exc})$, as shown in Fig. 1(b)
(more information is provided in the Supplemental Material).
We assume the two superconductors ($L$ and $R$) to have the same zero-
temperature energy gap ($\Delta_{L,0}=\Delta_{R,0}=\Delta_{0}$). Therefore, no
thermo-electric effect occurs for $h_{\rm exc}=0$ (Marchegiani  _et al._ ,
2020a, b).
Charge thermoelectric effect. Here, we investigate thermoelectric effects in
the $S-I-S_{m}$ junction as a function of thermal bias and exchange field.
Typical current-voltage characteristics for different values of $h_{\rm exc}$
at $T_{L}=0.7T_{c}$ and $T_{R}=0.01T_{c}$ are shown in Fig. 2(a). For $h_{\rm
exc}=0$ (black solid line), the system is dissipative [$I_{q}(V)V>0$], and the
current displays sub-gap peaks at $V=\pm
V_{p}=\pm|\Delta(T_{L})-\Delta(T_{R})|/|e|$. It is enough a weak exchange
field ($h_{\rm exc}=0.1\Delta_{0}$, magenta solid line) to observe the
splitting of the peaks at $V=\pm
V_{p}^{\mp}=\pm|\Delta(T_{L})-\Delta(T_{R})\pm h_{\rm exc}|/|e|$. For larger
values of the exchange field ($h_{exc}\gtrsim 0.2\Delta_{0}$), a
thermoelectric power is generated [$I_{q}(V)V<0$] 555Thermoelectricity is
identified by the presence of IV in the top-left or bottom-right quadrant of
Fig. 2(a), where necessary $I(V)V<0$. The exchange field in the presence of a
large thermal bias generates spontaneously a thermo-spin effect (discussed in
more details below). The spin-splitting in $S_{m}$ is fundamental to activate
this effect when $\Delta_{L,0}=\Delta_{R,0}$. Indeed, the spin-splitting in
the right electrode DoS reduces effectively the gap, and, at the same time,
localizes purely spin-polarized states between $\Delta_{R}-h_{exc}$ and
$\Delta_{R}+h_{exc}$. These combined mechanisms determine a spontaneous spin-
Seebeck effect, leading to thermoelectricity. When thermoelectricity is
generated, we can observe Absolute Negative Conductance (ANC), that is
$G(V)=I_{q}(V)/V<0$. In the limit of $V\to 0$, the conductance can be
approximated by

$\frac{G_{0}}{G_{T}}\simeq-\Delta_{0}^{2}\sum_{\sigma}\\!\int_{0}^{\infty}\\!\\!\\!\\!\\!\\!\\!dE\frac{N_{L}(E)f_{R}(E,T_{R})}{[(E+i\Gamma+\sigma
h_{\rm exc})^{2}-\Delta_{0}^{2}]^{3/2}}.$ (4)

See for instance the dashed violet line ($I_{q}=G_{0}V$) for
$h_{exc}=0.2\Delta_{0}$ in Fig. 2(a). We note that the above expression is
valid only for $\Delta(T_{L})>\Delta(T_{R})-h_{\rm exc}$, that is when
thermoelectricity appears. Figure 2(b) displays $G_{0}$ computed through
numerical differentiation of $I_{q}$ as a function of $T_{L}/T_{c}$ and
$h_{\rm exc}/\Delta_{0}$. We can distinguish thermoelectric ($G_{0}<0$, blue)
and dissipative ($G_{0}>0$, red) areas. For a fixed $h_{\rm exc}$,
thermoelectricity is present only in a limited range of $T_{L}$. In
particular, the maximum value of $T_{L}$ providing a thermoelectric effect is
due to the condition $\Delta(T_{L})>\Delta(T_{R})-h_{\rm exc}$ (dashed white
line), which corresponds to the requirement to have the hot electrode with the
largest ”effective” gap (Marchegiani  _et al._ , 2020a, b).
We now consider thermoelectricity at internal peak voltages
$V=\pm{V_{p}^{-}}=\pm|\Delta(T_{L})-(\Delta(T_{R})-h_{\rm exc})|/|e|$, where
typically the thermopower is maximal. Being $\Delta(T_{L})<\Delta(T_{R})$, we
note that $|\rm{V_{p}^{-}}|$ can be increased by raising $h_{\rm exc}$. As a
consequence, the generated maximum thermopower $-I_{q}(V_{p}^{-})V_{p}^{-}$ is
expected to increase. Note that $h_{\rm exc}$ cannot be freely increased, due
to the Chandrasekhar-Clogston limit ($h_{\rm exc}\leq\Delta_{0}/\sqrt{2}$)
(Clogston, 1962; Chandrasekhar, 1962; Maki, 1964).
A typical thermoelectric figure of merit is the Seebeck voltage ($V_{S}$),
which represents the bias that stops the thermocurrent [see the orange circles
in Fig. 2(a)]. This voltage grows monotonically with $h_{\rm exc}$, as shown
in Fig. 2(c) for different values of $T_{L}$, and its maximum is limited again
by the Chandrasekhar-Clogston limit (See the Supplemental Material). Moreover,
$|V_{S}|$ increases by lowering $T_{L}$, which corresponds to a _decrease_ of
the thermal gradient. This odd behavior shows that thermoelectricity is
directly connected to the nonlinearity of the effect (Marchegiani  _et al._ ,
2020a), which is different from thermoelectricity in linear regime Ozaeta  _et
al._ (2014). Another significant figure of merit is the ratio between the
maximum thermoelectric current and the corresponding voltage, i.e., $G^{\rm
max}\sim I^{\pm}(V_{p}^{-})/V_{p}^{-}$ [see aquamarine dotted line for $h_{\rm
exc}=0.4\Delta_{0}$ in Fig. 2(a)]. This quantity plays a crucial role when the
thermoelectric element is connected to a load. More precisely, $-G^{\rm max}$
represents the maximal conductance of the load supported by the thermoelectric
junction, such as no net thermopower can be generated if the load conductance
is bigger than $-G^{\rm max}$ Marchegiani  _et al._ (2020b). Figure 2(d) shows
$G^{\rm max}$ as a function of $h_{\rm exc}$ and $T_{L}$. We can identify
again thermoelectric ($G^{\rm max}\\!<0$, blue) and dissipative ($G^{\rm
max}>0$, red) regimes. We find that the temperature range where the junction
is thermo-active widens by increasing $h_{\rm exc}$. Indeed, for a given
$h_{\rm exc}$, the maximum value of $T_{L}$ is still limited by the above-
mentioned relation for $G_{0}$ [dashed white line of Fig.2(b)-2(d)]. For low
temperatures, thermoelectricity disappears arising from the nonlinear nature
of the effect in temperature 666The lower limit in temperature is nonuniversal
being related to the chosen $\Gamma$ parameter.. Some other differences
between Fig. 2(b) and Fig. 2(d) are discussed in more details in the
Supplemental Material.

Figure 2: (a) Quasiparticle current-voltage ($I_{q}V$) characteristics are
shown for different values of $h_{\rm exc}/\Delta_{0}$. In the linear regime,
the IV curves can be approximated by the zero-bias conductance ($G_{0}$) of
Eq. (4) (dashed violet line), while, in the nonlinear regime in bias, the
$G^{max}$ evaluated at the thermoelectric peaks (dotted aquamarine line). The
Seebeck voltages ($\pm V_{S}$) (orange circles) and the matching peaks values
($\pm V_{p}$) are reported for $h_{ex}=0.4\Delta_{0}$. (b) In the linear
regime, the zero-bias conductance ($G_{0}$) as a function of $h_{\rm exc}$ and
$T_{L}$ is shown, distinguishing thermo-active areas (blue tones) from
dissipative ones (red tones). The negative conductance is limited in
temperature by the relation $\Delta(T_{L})>\Delta(T_{R})-h_{\rm exc}$ (dashed
white line) (c) Absolute value of the Seebeck voltage as a function of $h_{\rm
exc}$ is shown for different values of the thermal bias ($T_{L}$). (d) Maximum
conductance evaluated at the matching peak ${V_{p}}^{-}$ is displayed as a
function of $T_{L}$ and $h_{\rm exc}$ in the nonlinear regime. The blue tones
are linked to the negative conductance (thermoelectricity), while the red area
is referred to the positive one (dissipation).

Spin thermoelectric effect. Here, we investigate the current spin-
polarization, which represent another peculiar consequence of spin symmetry
breaking in the setup considered.

Figure 3: (a) The charge current ($I_{q}$) is separated into up
($I_{\uparrow}$) and down ($I_{\downarrow}$) spin components for $h_{\rm
exc}=0.4\Delta$. Only one spin component ($I_{\downarrow}$ for $V>0$,
$I_{\uparrow}$ for $V<0$) generates a thermo-active peak. The spin current
($I_{s}$) is reported in the dashed orange line. (b) The IV characteristics of
the charge current ($I_{q}$) and spin current ($I_{s}$) and the curve of the
spin current generation efficiency ($SGE$) are displayed for $h_{\rm
exc}=0.4\Delta$. The yellow area between the Seebeck voltages ($\pm V_{s}$)
highlights the thermo-active spin current generation, while in the violet area
the system is dissipative behaving as a spin-filter. (c) The spin current
evaluated at $+{{V}_{p}}^{-}$ is reported as a function of $T_{L}$ and $h_{\rm
exc}$. Fixing $T_{L}$ (colored lines), the maximal value of $I_{s}$ is reached
for the optimal value of $h_{\rm exc}$ such as $\Delta_{L}=\Delta_{R}-h_{\rm
exc}$ (dashed white line). With a black dashed line, we show the approximation
of $h^{t}_{\rm exc}$ calculated in the Supplemental Material. The violet
transparent zone highlights the dissipative behavior of the system. (d) The
maximum spin current, evaluated at the matching peak voltage
($+{{V}_{p}}^{-}$), corresponding to the colored cuts of the previous panel(c)
(dashed lines for the dissipative regime and solid lines for the
thermoelectric regime), is shown as a function of $h_{\rm exc}$ for different
values of $T_{L}$.

Figure 3(a) shows typical charge current ($I_{q}$, aquamarine line) and spin
current ($I_{s}$, orange dashed line) as a function of $V$ in the presence of
thermoelectricity. In addition, we consider the two spin components of the
current, i.e., $I_{\uparrow}$ (spin up, blue line) and $I_{\downarrow}$ (spin
down, light green line). At $V=0$, $I_{q}(V=0)=0$, due to PH symmetry. Still,
the two spin components are finite, and exactly opposite, since $h_{\rm
exc}\neq 0$ [and so $N_{R_{\uparrow}}(E)\neq N_{R_{\uparrow}}(-E)$].
Therefore, one of the spin components is thermoactive [in Fig. 3(a)
$I_{\uparrow}(V)V<0$ for $V<0$, $I_{\downarrow}(V)V<0$ for $V>0$], while the
other one is dissipative. In this respect, the system produces spin-polarized
thermoelectricity when the thermoactive component is larger than the
dissipative one. Interestingly, the Seebeck voltage represents the bias where
the thermoactive and dissipative components compensate, i.e.,
$I_{\uparrow}(+V_{S})=-I_{\downarrow}(+V_{S})$ such as $I_{q}(\pm V_{S})=0$,
thus obtaining a purely spin current. Moreover, the maximum spin current is
located at the internal matching peaks ($\pm\rm{{V_{p}}^{-}}$).
It is convenient to introduce a new figure of merit, which compares the spin
current generated by the system with respect to all the tunneling carriers
moving in the junction. Indeed, the Spin current Generation Efficiency is
defined as $SGE=I_{S}/(|I_{\uparrow}|+|I_{\downarrow}|)$. Note that when the
two spin components flow in opposite directions, $|I_{s}|>|I_{q}|$ and the
system is globally thermo-active [$I_{q}(V)V<0$], necessarily $SGE=1$. Figure
3(b) displays $SGE$ as a function of $V$, for the physical situations of Fig.
3(a). For $|V|>V_{S}$, where the junctions are dissipative (violet area), we
can distinguish two behaviors. The first is characterized by $SGE=1$ when a
thermoactive spin component but lower than the dissipative one is still
present and hence no charge thermopower is generated. In the second case,
$SGE<1$. Here, both spin components are dissipative and the total charge
current flows in the same direction of the bias. The current is still spin-
polarized, but $|I_{s}|\leq|I_{q}|$ similarly to a spin-polarizer (Ozaeta  _et
al._ , 2014; Linder and Bathen, 2016; Bergeret  _et al._ , 2018). For large
$|V|$, the total current is almost independent of the spin-splitting, and the
two components give the same contribution, leading to $I_{s}\to 0$.
In Fig. 3(c), we analyse the maximal spin current evaluated at the peak
$I_{s}(\pm V_{p}^{-})$ as a function of $h_{\rm exc}$ and $T_{L}$. At fixed
$h_{\rm exc}$, $I_{s}$ increases with the thermal bias $\delta
T=T_{L}-T_{R}\sim T_{L}$ 777However, to remain thermoactive we need to stay
under the white dashed line.. By contrast, for a given thermal gradient
($\delta T$, colored cuts), the spin current is non-monotonic in $h_{\rm
exc}$, reaching the maximum value for $h_{\rm
exc}=\Delta(T_{R})-\Delta(T_{L})$. This condition coincides with the threshold
value of the exchange field [$h_{\rm exc}^{t}(T_{L})$] to generate
thermoelectricity. The threshold value can be estimated as $h^{t}_{\rm
exc}=\Delta\left(1-\tanh\left[1.74\sqrt{(T_{c}/T_{L})-1}\right]\right)$ (see
Supplemental Material). Its approximate value is displayed in Fig. 3(c) with a
black dashed line. In particular, for a fixed thermal bias $T_{L}$,
thermoelectricity occurs only for $h_{\rm exc}\geq h_{\rm exc}^{t}$ (solid
part of the colored cuts), otherwise, it is dissipative (dashed part).
Finally, Fig. 3(d) shows $I_{s}(\pm V_{p})$ as a function of $h_{exc}$ for
selected values of $T_{L}$. In particular, $I_{s}(V_{p}^{-})$ grows with
$T_{L}$, reaching its maximum at $h^{t}_{\rm exc}$. We observe that the
maximum value of the spin current is of the order of $I^{max}_{s}\sim
0.5G_{T}\Delta_{0}/e$. Assuming the aluminium gap of $\Delta_{0}=200~{}\mu eV$
and tunneling conductance of $G_{T}=0.1~{}mS$, we expect spin currents of the
order of 10 nA.
Conclusions. In summary, we discussed the nonlinear thermo-spin effect,
generated in a thermally-biased $S-I-S_{m}$ Josephson junction in the presence
of an exchange interaction. This effect is generated by the spontaneous
particle-hole symmetry breaking activated by spin-splitting in $S_{m}$. By
exploiting the two spin current components that the system naturally drives in
opposite direction, we observe the coexistence of one spin thermo-active with
one dissipative component. Notably, the spin thermoelectricity is relevant
when the thermoactivated component dominates, thus obtaining
$|I_{s}|>|I_{q}|$. Our results suggest interesting applications in
thermoelectricity (Benenti  _et al._ , 2017), low-dissipative and thermoactive
spintronics (Žutić  _et al._ , 2004; Linder and Bathen, 2016), and radiation
detection (Geng  _et al._ , 2020).

We acknowledge the European Research Council under Grant Agreement No.
899315-TERASEC, and the EU’s Horizon 2020 research and innovation program
under Grant Agreement No. 800923 (SUPERTED) and No. 964398 (SUPERGATE) for
partial financial support. A. B. acknowledges the CNR-CONICET cooperation
program Energy conversion in quantum nanoscale hybrid devices, the Royal
Society through the International Exchanges between the UK and Italy (No. IEC
R2 192166), and the SNS-WIS joint lab QUANTRA funded by the Italian Ministry
of Foreign Affairs and International Cooperation.

## References

  * Žutić  _et al._ (2004) I. Žutić, J. Fabian, and S. Das Sarma, Rev. Mod. Phys. 76, 323 (2004).
  * Linder and Robinson (2015) J. Linder and J. W. A. Robinson, Nat. Phys. 11, 307 (2015).
  * Tedrow and Meservey (1971) P. M. Tedrow and R. Meservey, Phys. Rev. Lett. 27, 919 (1971).
  * Meservey  _et al._ (1975) R. Meservey, P. M. Tedrow, and R. C. Bruno, Phys. Rev. B 11, 4224 (1975).
  * Takahashi  _et al._ (1999) S. Takahashi, H. Imamura, and S. Maekawa, Phys. Rev. Lett. 82, 3911 (1999).
  * Jedema  _et al._ (1999) F. J. Jedema, B. J. van Wees, B. H. Hoving, A. T. Filip, and T. M. Klapwijk, Phys. Rev. B 60, 16549 (1999).
  * Takahashi  _et al._ (2001) S. Takahashi, T. Yamashita, T. Koyama, S. Maekawa, and H. Imamura, J. Appl. Phys. 89, 7505 (2001).
  * Yamashita  _et al._ (2002) T. Yamashita, S. Takahashi, H. Imamura, and S. Maekawa, Phys. Rev. B 65, 172509 (2002).
  * Maekawa  _et al._ (2002) S. Maekawa, S. Takahashi, and H. Imamura, J. Phys. D 35, 2452 (2002).
  * Giazotto and Taddei (2008) F. Giazotto and F. Taddei, Phys. Rev. B 77, 132501 (2008).
  * Hao  _et al._ (1990) X. Hao, J. S. Moodera, and R. Meservey, Phys. Rev. B 42, 8235 (1990).
  * Giazotto (2009) F. Giazotto, Appl. Phys. Lett. 95, 042503 (2009).
  * Bergeret  _et al._ (2012) F. S. Bergeret, A. Verso, and A. F. Volkov, Phys. Rev. B 86, 214516 (2012).
  * De Simoni  _et al._ (2018) G. De Simoni, E. Strambini, J. S. Moodera, F. S. Bergeret, and F. Giazotto, Nano Lett. 18, 6369 (2018).
  * Virtanen and Heikkilä (2004) P. Virtanen and T. T. Heikkilä, J. Low Temp. Phys. 136, 401 (2004).
  * Titov (2008) M. Titov, Phys. Rev. B 78, 224521 (2008).
  * Jacquod and Whitney (2010) P. Jacquod and R. S. Whitney, EPL 91, 67009 (2010).
  * Kalenkov and Zaikin (2017) M. S. Kalenkov and A. D. Zaikin, Phys. Rev. B 95, 024518 (2017).
  * Heidrich and Beckmann (2019) J. Heidrich and D. Beckmann, Phys. Rev. B 100, 134501 (2019).
  * Eom  _et al._ (1998) J. Eom, C.-J. Chien, and V. Chandrasekhar, Phys. Rev. Lett. 81, 437 (1998).
  * Parsons  _et al._ (2003) A. Parsons, I. A. Sosnin, and V. T. Petrashov, Phys. Rev. B 67, 140502 (2003).
  * Jiang and Chandrasekhar (2005) Z. Jiang and V. Chandrasekhar, Phys. Rev. B 72, 020502 (2005).
  * Giazotto  _et al._ (2015) F. Giazotto, P. Solinas, A. Braggio, and F. S. Bergeret, Phys. Rev. Applied 4, 044016 (2015).
  * Tan  _et al._ (2021) Z. B. Tan, A. Laitinen, N. S. Kirsanov, P. J. Hakonen, A. Galda, V. M. Vinokur, M. Haque, A. Savin, D. S. Golubev, G. B. Lesovik, and P. J. Hakonen, Nat. Comm. 12, 138 (2021).
  * Hussein  _et al._ (2019) R. Hussein, M. Governale, S. Kohler, W. Belzig, F. Giazotto, and A. Braggio, Phys. Rev. B 99, 075429 (2019).
  * Kirsanov  _et al._ (2019) N. S. Kirsanov, Z. B. Tan, D. S. Golubev, P. J. Hakonen, and G. B. Lesovik, Phys. Rev. B 99, 115127 (2019).
  * Blasi  _et al._ (2020a) G. Blasi, F. Taddei, L. Arrachea, M. Carrega, and A. Braggio, Phys. Rev. Lett. 124, 227701 (2020a).
  * Blasi  _et al._ (2020b) G. Blasi, F. Taddei, L. Arrachea, M. Carrega, and A. Braggio, Phys. Rev. B 102, 241302 (2020b).
  * Heikkilä  _et al._ (2018) T. T. Heikkilä, R. Ojajärvi, I. J. Maasilta, E. Strambini, F. Giazotto, and F. S. Bergeret, Phys. Rev. Applied 10, 034053 (2018).
  * Machon  _et al._ (2013) P. Machon, M. Eschrig, and W. Belzig, Phys. Rev. Lett. 110, 047002 (2013).
  * Ozaeta  _et al._ (2014) A. Ozaeta, P. Virtanen, F. S. Bergeret, and T. T. Heikkilä, Phys. Rev. Lett. 112, 057001 (2014).
  * Linder and Bathen (2016) J. Linder and M. E. Bathen, Phys. Rev. B 93, 224509 (2016).
  * Bergeret  _et al._ (2018) F. S. Bergeret, M. Silaev, P. Virtanen, and T. T. Heikkilä, Rev. Mod. Phys. 90, 041001 (2018).
  * Kolenda  _et al._ (2017) S. Kolenda, C. Sürgers, G. Fischer, and D. Beckmann, Phys. Rev. B 95, 224505 (2017).
  * Shapiro  _et al._ (2017) D. S. Shapiro, D. E. Feldman, A. D. Mirlin, and A. Shnirman, Phys. Rev. B 95, 195425 (2017).
  * Sánchez and López (2016) D. Sánchez and R. López, Compendus R. Phys. 17, 1060 (2016).
  * Benenti  _et al._ (2017) G. Benenti, G. Casati, K. Saito, and R. S. Whitney, Phys. Rep. 694, 1 (2017).
  * Marchegiani  _et al._ (2020a) G. Marchegiani, A. Braggio, and F. Giazotto, Phys. Rev. Lett. 124, 106801 (2020a).
  * Marchegiani  _et al._ (2020b) G. Marchegiani, A. Braggio, and F. Giazotto, Phys. Rev. B 101, 214509 (2020b).
  * Tokuyasu  _et al._ (1988) T. Tokuyasu, J. A. Sauls, and D. Rainer, Phys. Rev. B 38, 8823 (1988).
  * Hao  _et al._ (1991) X. Hao, J. S. Moodera, and R. Meservey, Phys. Rev. Lett. 67, 1342 (1991).
  * Strambini  _et al._ (2017) E. Strambini, V. N. Golovach, G. De Simoni, J. S. Moodera, F. S. Bergeret, and F. Giazotto, Phys. Rev. Materials 1, 054402 (2017).
  * Tedrow  _et al._ (1986) P. M. Tedrow, J. E. Tkaczyk, and A. Kumar, Phys. Rev. Lett. 56, 1746 (1986).
  * Moodera  _et al._ (1988) J. S. Moodera, X. Hao, G. A. Gibson, and R. Meservey, Phys. Rev. Lett. 61, 637 (1988).
  * Tinkham (2004) M. Tinkham,  _Introduction to superconductivity_ (Dover Publications, Mineola, New York, 2004).
  * Note (1) We assume the insulating barrier is sufficiently opaque suppressing the Josephson coupling or the in-plane magnetic field suppresses it by Franhoufer interference Marchegiani  _et al._ (2020c); Barone (1992).
  * Shapiro  _et al._ (1962) S. Shapiro, P. H. Smith, J. Nicol, J. L. Miles, and P. F. Strong, IBM J. Res. Dev. 6, 34 (1962).
  * Note (2) For simplicity, we assume that the chemical potential for the two spin-species is equal $\mu_{\uparrow}=\mu_{\downarrow}$ (no spin accumulation). This approximation is valid for sufficiently opaque tunneling barrier.
  * Note (3) The self-consistency relation is Giazotto (2009); Giazotto  _et al._ (2015):  $\mathop{log}\nolimits\left(\frac{\Delta_{0}}{\Delta_{\alpha}}\right)=\intop\ilimits@_{0}^{\hbar\omega_{D}}dE\tmspace+{.1667em}\frac{f_{+}(E)-f_{-}(E)}{\sqrt{E^{2}+\Delta_{\alpha}^{2}}},$ (5) where $\Delta_{0}$ is the zero-temperature pairing potential in absence of the exchange field, $\omega_{D}$ is the Debye frequency and $f_{\pm}(E)=(1+\mathop{exp}\nolimits[(\sqrt{E^{2}+\Delta^{2}}\mp h_{\rm exc})/k_{B}T_{R}])^{-1}$ with $k_{B}$ Boltzmann constant.
  * Dynes  _et al._ (1984) R. C. Dynes, J. P. Garno, G. B. Hertel, and T. P. Orlando, Phys. Rev. Lett. 53, 2437 (1984).
  * Note (4) $G_{T}$ is defined as $G_{T}=G_{\uparrow}+G_{\downarrow}$. We assume that the two spin contributes are equal, such as $G_{\uparrow}=G_{\downarrow}$.
  * Note (5) Thermoelectricity is identified by the presence of IV in the top-left or bottom-right quadrant of Fig. 2(a), where necessary $I(V)V<0$.
  * Clogston (1962) A. M. Clogston, Phys. Rev. Lett. 9, 266 (1962).
  * Chandrasekhar (1962) B. S. Chandrasekhar, Appl. Phys. Lett. 1, 7 (1962).
  * Maki (1964) K. Maki, Prog. Theor. Phys. 31, 731 (1964).
  * Note (6) The lower limit in temperature is nonuniversal being related to the chosen $\Gamma$ parameter.
  * Note (7) However, to remain thermoactive we need to stay under the white dashed line.
  * Geng  _et al._ (2020) Z. Geng, A. P. Helenius, T. T. Heikkilä, and I. J. Maasilta, J. Low Temp. Phys. 199, 585 (2020).
  * Marchegiani  _et al._ (2020c) G. Marchegiani, A. Braggio, and F. Giazotto, Phys. Rev. Research 2, 043091 (2020c).
  * Barone (1992) A. Barone,  _Principles and Applications of Superconducting Quantum Interference Devices_ (Word Scientific, Singapore, 1992).

## Supplementary Information

##  I SYMMETRIES OF THE CHARGE CURRENT AND THE SPIN CURRENT

We discuss the symmetries of the charge current ($I_{q}$) and the spin current
($I_{s}$), which directly derive from the particle-hole (PH) symmetry of the
superconducting density of states (DoS). As discussed in the main text, the
spin-$\sigma$ (with $\sigma=\pm$) DoS of the $\alpha$ superconductor (with
$\alpha=L,R$) in the junction is Tinkham (2004)

$N_{\alpha\sigma}(E_{\alpha})=\left|\Re\\!\\!\left[\frac{E_{\alpha}+\sigma
h_{\alpha}+i\Gamma}{\sqrt{(E_{\alpha}+\sigma
h_{\alpha}+i\Gamma)^{2}-{\Delta_{\alpha}}^{2}}}\right]\\!\\!\right|,$ (S1)

where $h_{\alpha}$ is the exchange field, $E_{\alpha}=E-\mu_{\alpha}$ is the
quasiparticle energy measured with respect to the chemical potential
$\mu_{\alpha}$, and $\Delta_{\alpha}$ is the self-consistent gap for the
$\alpha$-lead, $\Gamma$ is the Dynes parameter. The full DoS is given by the
sum $N_{\alpha}(E)=\sum_{\sigma=\pm}N_{\alpha\sigma}(E)/2$. In the presence of
an exchange (or Zeeman) field $h_{\rm exc}$, even if the spin degeneracy of
$S_{m}$ is broken, the superconducting state of $S_{m}$ still satisfies the PH
symmetry, which implies the symmetry
$N_{\alpha\sigma}(E_{\alpha})=N_{\alpha\bar{\sigma}}(-E_{\alpha})$, where
$\bar{\sigma}=-\sigma$. Anyway, the PH symmetry of the BCS superconductors is
reflected in the BCS DoS by its even symmetry with respect to the energy, i.e.
$N_{\alpha}(E_{\alpha})=N_{\alpha}(-E_{\alpha})$.
In the following, we assume to have the exchange field _only_ in the right
lead, i. e., $h_{L}=0$, as discussed in the setup S-I-$S_{m}$, shown in Fig.1
of the main text. Correspondingly, the current for the $\sigma$-spin component
becomes

$I_{\sigma}(V)=\frac{G_{T}}{2e}\int_{-\infty}^{\infty}\\!dE~{}N_{L\sigma}(E-eV)N_{R\sigma}(E,h_{exc}){F}_{LR}(E,V),$
(S2)

as reported in the main text. In the above expression,
${F}_{LR}(E,V)=f_{L}(E-eV,T_{L})-f_{R}(E,T_{R})$ is the difference between the
Fermi-Dirac quasiparticle distributions of two electrodes
$f_{\alpha}(E,T_{\alpha})=[1+\exp(E/k_{B}T_{\alpha})]^{-1}$. Note that the PH
symmetry of the Fermi function with respect to the energy gives
$f_{\alpha}(E)=1-f_{\alpha}(-E)$ and determines also the symmetry
${F}_{LR}(-E,V)=-{F}_{LR}(E,-V)$. By changing the integration variable
$E\to-E^{\prime}$ in Eq. (S2) and using the symmetries of
$N_{\alpha\sigma}(E)$ and ${F}_{LR}(-E,V)$, we obtain

$I_{\sigma}(V)=-\frac{G_{T}}{2e}\int_{-\infty}^{\infty}\\!dE^{\prime}~{}N_{L\sigma}(E^{\prime}+eV)N_{R\bar{\sigma}}(E^{\prime},h_{exc}){F}_{LR}(E^{\prime},-V)=-I_{\bar{\sigma}}(-V),$
(S3)

where the second identity is obtained by direct comparison with Eq. (S2). By
exploiting the symmetry of $I_{\sigma}$, it is easy to show that charge
($I_{q}=\sum_{\sigma=\pm}I_{\sigma}$) and spin ($I_{s}=\sum_{\sigma=\pm}\sigma
I_{\sigma}$) currents are odd and even functions in $V$, respectively. For
instance, in the case of the symmetry with respect to the voltage bias $V$,
one gets

$\begin{split}I_{q}(-V)&=-\sum_{\sigma=\pm}I_{\bar{\sigma}}(V)=-I_{q}(V)\\\
I_{s}(-V)&=-\sum_{\sigma=\pm}\sigma
I_{\bar{\sigma}}(V)=\sum_{\sigma=\pm}\bar{\sigma}I_{\bar{\sigma}}(V)=I_{s}(V)\end{split}$
(S4)

Similarly, one can easily show that charge (spin) current is even (odd) with
respect to the exchange field $I_{q}(V,h_{exc})=I_{q}(V,-h_{exc})$
[$I_{s}(V,h_{exc})=-I_{s}(V,-h_{exc})$].

##  II Evaluation of $V_{p}^{-}$ and the $h_{exc}$ thermoelectrical threshold

In the presence of thermoelectricity, and assuming that the critical
temperature $T_{c}$ of the superconductors is $T_{c}\gtrsim T_{L}\gg T_{R}$,
the value of the voltage bias where the thermoelectric current is maximum
reads

$eV_{p}^{-}=\Delta(T_{L})-\left[\Delta(T_{R},h_{exc})-h_{exc}\right].$ (S5)

For low temperatures ($T_{R}\ll T_{c}$) and low exchange fields, we can
approximate $\Delta(T_{R},h_{exc})\approx\Delta_{0}$. We can estimate the
maximum value of the peak voltage ($V_{p}^{-,max}$) by considering the maximum
exchange field allowed, i.e., $h_{exc}\leq h_{exc}^{max}=\Delta_{0}/\sqrt{2}$,
given by the Chandrasekhar-Clogston limit (Clogston, 1962; Chandrasekhar,
1962). We obtain

$eV_{p}^{-,max}\approx\Delta_{0}(\tanh[1.74\sqrt{T_{c}/T_{L}-1}]-0.29),$ (S6)

where we approximate $\Delta(T_{L})=\Delta_{0}\tanh[1.74\sqrt{T_{c}/T_{L}-1}]$
Tinkham (2004), which overall differs only up to a few percent from the
results obtained through self-consistent calculation of $\Delta_{L}$. This
limit provides also a rough estimation of the Seebeck voltage values, being
$V_{S}\gtrsim V_{p}^{-,max}$. By using a similar approximation, we can compute
the minimum threshold value for $h_{\rm exc}$ necessary to generate
thermoelectricity in the system. We observe that the matching peak is located
at zero bias, i.e., $V_{p}^{-}=0$, when the system switches from dissipative
to thermoactive behavior. Thus, using again the Eq. (S5) we find

$h_{\rm exc}^{t}=\Delta_{R}(T_{R},h_{\rm exc}^{t})-\Delta_{L}(T_{L},h_{\rm
exc}=0).$ (S7)

Using previous approximations, we find

$\Delta_{0}(\tanh[1.74\sqrt{T_{c}/T_{L}-1}]-[1-h_{\rm exc}^{t}/\Delta_{0}])=0$
(S8)

and the threshold $h_{exc}^{t}$ becomes

$h_{\rm
exc}^{t}=\Delta_{0}\left[1-\tanh\left(1.74\sqrt{\frac{T_{C}}{T_{L}}-1}\right)\right].$
(S9)

The relation implies that by increasing the temperature ($T_{L}$), $h_{\rm
exc}^{t}$ rises as shown by the Seebeck voltage $V_{S}$ of Fig. 2(c) of the
main text. From the previous results, we expect that the Seebeck voltage
$V_{S}$ depends linearly on $h_{\rm exc}$ and $h_{\rm exc}^{t}(T_{L})$. We can
provide a very rough estimation of the Seebeck voltage with

$eV_{S}^{-}\approx h_{\rm exc}-h_{\rm exc}^{t}(T_{L})=h_{\rm
exc}-\Delta\left[1-\tanh\left(1.74\sqrt{\frac{T_{C}}{T_{L}}-1}\right)\right].$
(S10)

##  III The negative conductance in linear and non linear regime

In the main text, we discussed the negative conductance ($G_{0}$) in the
linear regime and the maximal conductance ($G^{\rm max}$) in the nonlinear
regime as a function of $h_{\rm exc}$ and $T_{L}$ [Fig. S1(a-b)]. We observed
the condition for thermoelectricity (blue tones) and its limits in temperature
at fixed $h_{\rm exc}$.
Here, we focus on the differences that appear between the linear and nonlinear
regimes. In particular, observing the central area of both figures
($T_{L}/T_{c}\sim 0.4$), the system is dissipative $G_{0}>0$ around $V\sim 0$
[Fig. S1(a)], while it exhibits thermoelectric effects $G^{\rm max}<0$ at the
matching peak $V=\pm V_{p}^{-}$ [Fig. S1(b)].
For clarity, we report IV characteristics of the charge current at fixed
temperature ($T_{L}/T_{C}=0.4$) for different values of $h_{\rm exc}$ [Fig.
S1(c)], as indicated by the colored dots in the Fig. S1(a). We can note the
dissipative behavior [$I_{q}(V)/V>0$] around $V\sim 0$ and thermoelectricity
[$I_{q}(V)/V<0$] around $\pm V_{p}^{-}$. This difference between linear and
nonlinear in bias behavior has been reported also in absence of $h_{\rm exc}$
for nonlinear thermoelectric systems Marchegiani  _et al._ (2020b). For a more
complete discussion, we show also the spin current for the same values. Note
that the curves display a different functional shape with respect to the case
analyzed in the main text at $T_{L}/T_{c}=0.7$ and for $h_{\rm
exc}=0.4\Delta_{0}$, which represents the typical behavior [Inset of Fig.
1(b)].
We observe that the reverse behavior is also observed. Indeed, for high values
of $h_{\rm exc}$ at $T_{L}/T_{c}\sim 0.2$ , the thermoelectricity is present
only in the linear regime ($V\sim 0$) [Fig. S1(a)]. In this case, the
thermoelectric currents are very small and strongly dependent on the Dynes
parameter representing a non typical feature.

Figure S1: (a) The negative conductance $G_{0}$ is reported as a function of
$T_{L}$ and $h_{\rm exc}$. Blue tones (red tones) show a thermoelectric
(dissipative) behavior. The colored points are referred to the IV
characteristics of Fig. S1(c) at $T_{L}/T_{c}=0.4$ and $h_{\rm exc}$ fixed.
(b) The maximal conductance $G^{\rm max}$ evaluated at $\pm V_{p}^{-}$ is
reported as a function of $T_{L}$ and $h_{\rm exc}$. In the same way, blue
tones (red tones) show thermoelectricity (dissipative regime). The colored
points are referred to the IV characteristics of Fig. S1(c) at $T_{L}$ and
$h_{\rm exc}$ fixed. (c) The $I_{q}(V)$ characteristics of the charge current
are shown at $T_{L}/T_{C}=0.4$ for different values of $h_{\rm exc}$. (d) The
$I_{s}(V)$ characteristics of the spin current are displayed at
$T_{L}/T_{C}=0.4$ for different values of $h_{\rm exc}$.