It turns out that one can get a stronger result than demanded in the question: compute $\Delta(a,b)$ for any $a,b \in S_{d-1}$, perpendicular or not. Indeed,

Claim. If $f(\rho)=a_0 + a_1 \rho + a_2 \rho^2 + a_3 \rho^3 + \mathcal O(\rho^4)$, then for every $u,v \in S_{d-1}$, and $x$ be uniform on the sphere, then we have the following approximation \begin{eqnarray} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] &= -\frac{a_2^2}{d^2}\\ &+(\dfrac{a_1^2}{d}+\frac{6a_1a_3}{d^2})u^\top v\\ & + \dfrac{2a_1^2}{d^2}(u^\top v)^2\\ &+\mathcal O(\dfrac{1}{d^3}). \end{split} \end{eqnarray} In paricular, if $u^\top v = 0$, then $$ \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] = -\frac{a_2^2}{d^2}+\mathcal O(\frac{1}{d^3}), $$ as proved by Iosif Pinelis (in the accepted answer).

For the proof, we will need the following result

Lemma. Let $u$ and $v$ be fixed and $x$ be uniformly random on the unit-sphere. Let $p$ and $q$ be nonnegative integers and define $c_{p,q}(u,v):=\mathbb E_x[(x^\top u)^p(x^\top v)^q]$. If $p$ and $q$ have different parities, then $c_{p,q}(u,v)=0$. Otherwise, we have the formula \begin{eqnarray} c_{p,q}(u,v) = \dfrac{p!q!\Gamma(\frac{d}{2})}{2^{p+q}\Gamma(\frac{d+p+q}{2})} \sum_t \dfrac{2^t}{t!(\frac{p-t}{2})!(\frac{q-t}{2})!}(u^\top v)^t, \end{eqnarray} where the sum is over all $t$ between $0$ and $p \land q$ inclusive, that have the same parity as $p$ and $q$. The formula is simplified in the table below for special values of $p$ and $q$. The above lemma is proved (ME link) here https://math.stackexchange.com/a/4004804/168758.

Proof of the claim. WLOG, let $a_0=0$. Thanks to the lemma, one may compute $$ (\mathbb E_x[f(x^\top u)])^2 = (\frac{a_2}{d}+\mathcal O(\frac{1}{d^2}))^2 = \frac{a_2^2}{d^2} + \mathcal O(\frac{1}{d^3}), $$ and similarly \begin{eqnarray*} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)] &= \sum_{p=0}^3\sum_{q=0}^3 a_p a_q\mathbb E_x[(x^\top u)^p(x^\top v)^q]+\mathcal O(\dfrac{1}{d^3})\\ &= \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d(d+2)}(u^\top v)^2 + \dfrac{6 a_1a_3}{d(d+2)}u^\top v+\mathcal O(\dfrac{1}{d^3}), %\\ &=a_0^2 + \frac{2a_0a_2}{d} + \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d^2}(u^\top v)^2 +\mathcal O(\dfrac{1}{d^4}), \end{split} \end{eqnarray*} and the claim follows after subtracting the previous display. $\quad\quad\quad\quad\Box$

It turns out that one can get a stronger result than demanded in the question: compute $\Delta(a,b)$ for any $a,b \in S_{d-1}$, perpendicular or not. Indeed,

Claim. If $f(\rho)=a_0 + a_1 \rho + a_2 \rho^2 + a_3 \rho^3 + \mathcal O(\rho^4)$, then for every $u,v \in S_{d-1}$, and $x$ be uniform on the sphere, then we have the following approximation \begin{eqnarray} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] &= -\frac{a_2^2}{d^2}\\ &+(\dfrac{a_1^2}{d}+\frac{6a_1a_3}{d^2})u^\top v\\ & + \dfrac{2a_1^2}{d^2}(u^\top v)^2\\ &+\mathcal O(\dfrac{1}{d^3}). \end{split} \end{eqnarray} In paricular, if $u^\top v = 0$, then $$ \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] = -\frac{a_2^2}{d^2}+\mathcal O(\frac{1}{d^3}), $$ as proved by Iosif Pinelis (in the accepted answer).

For the proof, we will need the following result

Lemma. Let $u$ and $v$ be fixed and $x$ be uniformly random on the unit-sphere. Let $p$ and $q$ be nonnegative integers and define $c_{p,q}(u,v):=\mathbb E_x[(x^\top u)^p(x^\top v)^q]$. If $p$ and $q$ have different parities, then $c_{p,q}(u,v)=0$. Otherwise, we have the formula \begin{eqnarray} c_{p,q}(u,v) = \dfrac{p!q!\Gamma(\frac{d}{2})}{2^{p+q}\Gamma(\frac{d+p+q}{2})} \sum_t \dfrac{2^t}{t!(\frac{p-t}{2})!(\frac{q-t}{2})!}(u^\top v)^t, \end{eqnarray} where the sum is over all $t$ between $0$ and $p \land q$ inclusive, that have the same parity as $p$ and $q$. The formula is simplified in the table below for special values of $p$ and $q$. The above lemma is proved (ME link) here https://math.stackexchange.com/a/4004804/168758.

Proof of the claim. WLOG, let $a_0=0$. Thanks to the lemma, one may compute $$ (\mathbb E_x[f(x^\top u)])^2 = (\frac{a_2}{d}+\mathcal O(\frac{1}{d^2}))^2 = \frac{a_2^2}{d^2} + \mathcal O(\frac{1}{d^3}), $$ and similarly \begin{eqnarray*} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)] &= \sum_{p=0}^3\sum_{q=0}^3 a_p a_q\mathbb E_x[(x^\top u)^p(x^\top v)^q]+\mathcal O(\dfrac{1}{d^3})\\ &= \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d(d+2)}(u^\top v)^2 + \dfrac{6 a_1a_3}{d(d+2)}u^\top v+\mathcal O(\dfrac{1}{d^3}) \\ &=(\dfrac{a_1^2}{d}+\frac{6a_1a_3}{d^2})u^\top v + \dfrac{2a_2^2}{d^2}(u^\top v)^2 +\mathcal O(\dfrac{1}{d^3}), \end{split} \end{eqnarray*} and the claim follows after subtracting the previous display. $\quad\quad\quad\quad\Box$

It turns out that one can get a stronger result than demanded in the question: compute $\Delta(a,b)$ for any $a,b \in S_{d-1}$, perpendicular or not. Indeed,

Claim. If $f(\rho)=a_0 + a_1 \rho + a_2 \rho^2 + a_3 \rho^3 + \mathcal O(\rho^4)$, then for every $u,v \in S_{d-1}$, and $x$ be uniform on the sphere, then we have the following approximation \begin{eqnarray} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] &= (\dfrac{a_1^2}{d}+\frac{6a_1a_3}{d^2})u^\top v\\ & + \dfrac{2a_1^2}{d(d+2)}(u^\top v)^2\\ &+\mathcal O(\dfrac{1}{d^3}). \end{split} \end{eqnarray} In paricular, if $u^\top v = 0$, then $$ \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] = \mathcal O(\frac{1}{d^3}), $$ as proved by Iosif Pinelis (in the accepted answer).

For the proof, we will need the following result

Lemma. Let $u$ and $v$ be fixed and $x$ be uniformly random on the unit-sphere. Let $p$ and $q$ be nonnegative integers and define $c_{p,q}(u,v):=\mathbb E_x[(x^\top u)^p(x^\top v)^q]$. If $p$ and $q$ have different parities, then $c_{p,q}(u,v)=0$. Otherwise, we have the formula \begin{eqnarray} c_{p,q}(u,v) = \dfrac{p!q!\Gamma(\frac{d}{2})}{2^{p+q}\Gamma(\frac{d+p+q}{2})} \sum_t \dfrac{2^t}{t!(\frac{p-t}{2})!(\frac{q-t}{2})!}(u^\top v)^t, \end{eqnarray} where the sum is over all $t$ between $0$ and $p \land q$ inclusive, that have the same parity as $p$ and $q$. The formula is simplified in the table below for special values of $p$ and $q$. The above lemma is proved (ME link) here https://math.stackexchange.com/a/4004804/168758.

Proof of the claim. WLOG, let $a_0=0$. Thanks to the lemma, one may compute $$ (\mathbb E_x[f(x^\top u)])^2 = (\frac{a_2}{d}+\mathcal O(\frac{1}{d^2}))^2 = \frac{a_2^2}{d^2} + \mathcal O(\frac{1}{d^3}), $$ and similarly \begin{eqnarray*} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)] &= \sum_{p=0}^3\sum_{q=0}^3 a_p a_q\mathbb E_x[(x^\top u)^p(x^\top v)^q]+\mathcal O(\dfrac{1}{d^3})\\ &= \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d(d+2)}(u^\top v)^2 + \dfrac{6 a_1a_3}{d(d+2)}u^\top v+\mathcal O(\dfrac{1}{d^3}), %\\ &=a_0^2 + \frac{2a_0a_2}{d} + \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d^2}(u^\top v)^2 +\mathcal O(\dfrac{1}{d^4}), \end{split} \end{eqnarray*} and the claim follows after subtracting the previous display. $\quad\quad\quad\quad\Box$

It turns out that one can get a stronger result than demanded in the question: compute $\Delta(a,b)$ for any $a,b \in S_{d-1}$, perpendicular or not. Indeed,

Claim. If $f(\rho)=a_0 + a_1 \rho + a_2 \rho^2 + a_3 \rho^3 + \mathcal O(\rho^4)$, then for every $u,v \in S_{d-1}$, and $x$ be uniform on the sphere, then we have the following approximation \begin{eqnarray} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] &= -\frac{a_2^2}{d^2}\\ &+(\dfrac{a_1^2}{d}+\frac{6a_1a_3}{d^2})u^\top v\\ & + \dfrac{2a_1^2}{d^2}(u^\top v)^2\\ &+\mathcal O(\dfrac{1}{d^3}). \end{split} \end{eqnarray} In paricular, if $u^\top v = 0$, then $$ \mathbb E_x[f(x^\top u)f(x^\top v)]-\mathbb E_x[f(x^\top u)]\mathbb E_x[f(x^\top v)] = -\frac{a_2^2}{d^2}+\mathcal O(\frac{1}{d^3}), $$ as proved by Iosif Pinelis (in the accepted answer).

For the proof, we will need the following result

Lemma. Let $u$ and $v$ be fixed and $x$ be uniformly random on the unit-sphere. Let $p$ and $q$ be nonnegative integers and define $c_{p,q}(u,v):=\mathbb E_x[(x^\top u)^p(x^\top v)^q]$. If $p$ and $q$ have different parities, then $c_{p,q}(u,v)=0$. Otherwise, we have the formula \begin{eqnarray} c_{p,q}(u,v) = \dfrac{p!q!\Gamma(\frac{d}{2})}{2^{p+q}\Gamma(\frac{d+p+q}{2})} \sum_t \dfrac{2^t}{t!(\frac{p-t}{2})!(\frac{q-t}{2})!}(u^\top v)^t, \end{eqnarray} where the sum is over all $t$ between $0$ and $p \land q$ inclusive, that have the same parity as $p$ and $q$. The formula is simplified in the table below for special values of $p$ and $q$. The above lemma is proved (ME link) here https://math.stackexchange.com/a/4004804/168758.

Proof of the claim. WLOG, let $a_0=0$. Thanks to the lemma, one may compute $$ (\mathbb E_x[f(x^\top u)])^2 = (\frac{a_2}{d}+\mathcal O(\frac{1}{d^2}))^2 = \frac{a_2^2}{d^2} + \mathcal O(\frac{1}{d^3}), $$ and similarly \begin{eqnarray*} \begin{split} \mathbb E_x[f(x^\top u)f(x^\top v)] &= \sum_{p=0}^3\sum_{q=0}^3 a_p a_q\mathbb E_x[(x^\top u)^p(x^\top v)^q]+\mathcal O(\dfrac{1}{d^3})\\ &= \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d(d+2)}(u^\top v)^2 + \dfrac{6 a_1a_3}{d(d+2)}u^\top v+\mathcal O(\dfrac{1}{d^3}), %\\ &=a_0^2 + \frac{2a_0a_2}{d} + \dfrac{a_1^2}{d}u^\top v + \dfrac{2a_2^2}{d^2}(u^\top v)^2 +\mathcal O(\dfrac{1}{d^4}), \end{split} \end{eqnarray*} and the claim follows after subtracting the previous display. $\quad\quad\quad\quad\Box$