Summary: Let \(f\) be a function on \(\mathbb{R}^2\) in the inhomogeneous Besov space \(\text{B}_{\infty ,1}^1(\mathbb{R}^2)\). For a pair \((A,B)\) of not necessarily bounded and not necessarily commuting self-adjoint operators, the function \(f(A,B)\) of \(A\) and \(B\) is introduced as a densely defined linear operator. It is shown that if \(1\le p\le 2\), \((A_1,B_1)\) and \((A_2,B_2)\) are pairs of not necessarily bounded and not necessarily commuting selfadjoint operators such that both \(A_1-A_2\) and \(B_1-B_2\) belong to the Schatten-von Neumann class \({\boldsymbol{S}}_p\) and \(f\in\text{B}_{\infty ,1}^1(\mathbb{R}^2)\), then the following Lipschitz type estimate holds: \[ \|f(A_1,B_1)-f(A_2,B_2)\|_{{\boldsymbol{S}}_p} \le \operatorname{const}\|f\|_{\text{B}_{\infty ,1}^1}\max \big \{\|A_1-A_2\|_{{\boldsymbol{S}}_p},\|B_1-B_2\|_{{\boldsymbol{S}}_p}\big \}. \]