There's a simple way by which "each round of DES algorithm is its own inverse". Consider round $n$ of DES as involving (almost only) a function $g_n$ with
$$g_n(L\mathbin\|R)=\bigl(L\oplus f(R,K_n)\bigr)\mathbin\|R$$
where $K_n$ is the 48-bit subkey for round $n$, function $f$ is the "cipher function" (given in the definition of DES), and $L$ and $R$ are 32-bit bitsrings forming a 64-bit block.
That function $g_n(L\mathbin\|R)$ verifies $g_n(g_n(L\mathbin\|R))=L\mathbin\|R$, as thought in the question; or in other words $g_n$ is an involution; or in yet other words $g_n\circ g_n$ is the identity function. Proof:
$$\begin{align}
g_n(g_n(L\mathbin\|R))&=g_n\Bigl(\bigl(L\oplus f(R,K_n)\bigr)\mathbin\|R\Bigr)\\
&=\Bigl(\bigl(L\oplus f(R,K_n)\bigr)\oplus f(R,K_n)\Bigr)\mathbin\|R\\
&=\Bigl(L\oplus\bigl(f(R,K_n)\oplus f(R,K_n)\bigr)\Bigr)\mathbin\|R\\
&=\left(L\oplus0^{32}\right)\mathbin\|R\\
&=L\mathbin\|R
\end{align}$$
That proof invokes the definition of $g_n$ (twice), associativity of $\oplus$, that $f$ is a function, that for all 32-bit $X$ it holds $X\oplus X=0^{32}$ (the bitstring of 32 zero bits), which is the neutral for $\oplus$.
DES encryption chains these 33 operations on 64-bit quantities:
$$\mathsf{IP}\,,\,g_1\,,\,\mathsf S\,,\,g_2\,,\,\mathsf S\,,\,\ldots\,,\,\mathsf S\,,\,g_{15}\,,\,\mathsf S\,,\,g_{16}\,,\,\mathsf{IP}^{-1}$$
where $\mathsf S$ is the "swap" involution defined by $\mathsf S(L\mathbin\|R)=R\mathbin\|L$, function $\mathsf{IP}$ is some permutation of bits, and $\mathsf{IP}^{-1}$ is the inverse permutation.
DES decryption chains these 33 operations on 64-bit quantities:
$$\mathsf{IP}\,,\,g_{16}\,,\,\mathsf S\,,\,g_{15}\,,\,\mathsf S\,,\,\ldots\,,\,\mathsf S\,,\,g_2\,,\,\mathsf S\,,\,g_1\,,\,\mathsf{IP}^{-1}$$
We see that DES encryption then decryption is the identity function: the $(34-j)^\text{th}$ operation of decryption cancels the $j^\text{th}$ operation of encryption:
- For $j=33$, because $\mathsf{IP}$ cancels $\mathsf{IP}^{-1}$.
- For $j=1$, because $\mathsf{IP}^{-1}$ cancels $\mathsf{IP}$.
- For other even $j$, because $g_{j/2}$ is an involution.
- For other (odd) $j$, because $\mathsf S$ is an involution.
Importantly, encryption and decryption use the very same structure, only the indexes (that, is, the order of the subkeys $K_n$) differ. That allow to use identical hardware or code for both encryption and decryption.
The usual definition of a round of a Feistel cipher includes the swap $\mathsf S$:
$$\begin{align}
g'_n(L\mathbin\|R)&=\mathsf S\bigl(g_n(L\mathbin\|R)\bigr)\\
&=R\mathbin\|\bigl(L\oplus f(R,K_n)\bigr)
\end{align}$$
and this function is not its own inverse per the sense in the question.
With that notation, DES encryption and decryption are the 18 operations
$$\mathsf{IP}\,,\,g'_1\,,\,g'_2\,,\,\ldots\,,\,g'_{15}\,,\,g_{16}\,,\,\mathsf{IP}^{-1}\\
\mathsf{IP}\,,\,g'_{16}\,,\,g'_{15}\,,\,\ldots\,,\,g'_2\,,\,g_1\,,\,\mathsf{IP}^{-1}$$
In that presentation, it is less apparent that decryption undoes encryption. But that still holds, and becomes obvious again if we expand the $g'_n$ into $g_n$ followed by $\mathsf S$.
When it is said that DES has 16 rounds, it is ignored $\mathsf{IP}$ and $\mathsf{IP}^{-1}$. And, for a strict Vulcan, only the 16th round (the one without swap) is its own inverse.
The DES specification, and many textbooks on Feistel ciphers, give indexes to $L$ and $R$ before and after round $n$ of encryption, with $L_0\mathbin\|R_0$ the plaintext after $\mathsf{IP}$. The above equation for $g'_n$ can then be written as:
$$\begin{align}
L_n&\gets R_{n-1}\\
R_n&\gets L_{n-1}\oplus f(R_{n-1},K_n)
\end{align}$$
DES (which as $m=16$ rounds) and some textbooks consider the ciphertext (before $\mathsf{IP}^{-1}$) to be $L_m\mathbin\|R_m$, and specialize the last encryption round's equations:
$$\begin{align}
L_m&\gets L_{m-1}\oplus f(R_{m-1},K_m)\\
R_m&\gets R_{m-1}
\end{align}$$
There are other conventions around: some texts use the same equations for all rounds. Then some add a final swap, or consider the ciphertext to be $R_m\mathbin\|L_m$; others consider $L_m\mathbin\|R_m$ as the ciphertext (that later kind does not obtain the same ciphertext as DES).
As to decryption, DES uses the same round equations for encryption and decryption except for the numbering of subkeys, that is for $1\le n<m$:
$$\begin{align}
L_n&\gets R_{n-1}&&&L_m&\gets L_{m-1}\oplus f(R_{m-1},K_1)\\
R_n&\gets L_{n-1}\oplus f(R_{n-1},K_{16-n})&&&R_m&\gets R_{m-1}
\end{align}$$
with the decryption's $L_0\mathbin\|R_0$ defined as the encryption's $L_n\mathbin\|R_n$, from which it follows that for $1\le n<m$, the decryption's $L_n\mathbin\|R_n$ is the encryption's $R_{m-n}\mathbin\|L_{m-n}$ (notice the inversion), and the decryption's $L_m\mathbin\|R_m$ is the encryption's $L_0\mathbin\|R_0$.
Other texts use the same naming for equal variables in decryption and encryption, and define different round equations for decryption. For regular round structure that gives:
$$\begin{align}
R_{n-1}&\gets L_n\\
L_{n-1}&\gets R_n\oplus f(L_n,K_n)
\end{align}$$
Footnote per comment: There is always some form of swap between rounds of a Feistel cipher, so that the fraction of the state that did not change in a round changes in the next round. That's essential for security (not for decryption to work). A round is its own inverse only if we exclude the swap as part of its definition.
