It is known that if $B=\prod_{i=1}^{\infty} \frac{z-z_i}{1-\overline{z_i}z}$ where $z_1, z_2,...$ are all distinct and $\sum_{i=1}^{\infty}(1-|z_n|)<\infty$, then $(BH^2)^{\perp}=\overline{span}\{\frac{1}{1-\overline{z_i}z} : i\in\mathbb{N}\}$
Now if we start with $f\in (BH^2)^{\perp}\bigcap H^2$, then $<Bz^n, f>=0$ $\forall n=0,1,...$
That means, $\frac{1}{2\pi}\int_0^{2\pi} B(e^{i\theta}) e^{in\theta} \overline{f(e^{i\theta})} d\theta =0$ $\forall n=0,1,...$
$\Rightarrow$ $B\overline{f}\in H^2$ $\Big[$ since, $B \overline{f} \in L^2$ and additionally, all the negative Fourier coefficients are zero. $\Big]$
Therefore, we can conclude $B h \in H^2$ where $h(z)= \frac{1}{1-z_i \overline{z}}$
But this is not true. For example, if we take $B(z)= \frac{z-\frac{1}{2}}{1-\frac{z}{2}}$, then $\frac{z-\frac{1}{2}}{1-\frac{z}{2}} \frac{1}{1-\frac{\overline{z}}{2}}$ is not analytic.
I am sure I have misunderstood something. I would appreciate it if you could assist me in any way.
This is very subtle!
You are right that $B \overline{f} \in H^2$ but, up to this point, you have treated $B\overline{f}$ as a boundary function. Therefore, you need to find an analytic $F$ on the unit disk with boundary value equal to $B\overline{f}$.
To illustrate this point, let us look at $h(z) = \frac{1}{1 - \frac{\overline{z}}{2}}$ as in your example. This is not analytic! The analytic function with the same boundary value as $Bh$ is, in fact, $F(z) = \frac{z - \frac{1}{2}}{1 - \frac{z}{2}} \frac{1}{1 - \frac{1}{2z}} = \frac{2z}{2-z}$.