Your function of p, call it F(p), if I understand you aright is F(p) = sum[all orientations and positions of polyominos X which contain the origin] p^|X| * q^|neighborset(X)| where |S| denotes the cardinality of [number of grid squares in] set S and where q=1-p. (Special case: when X is the emptyset then neighborset(X) also is empty.) It likely is best to think about the double generating function F(p,q) i.e. not assuming q+p=1 until afterwards. It is claimed (wikipedia) that the threshhold for square lattice site percolation is pcrit=0.5927460 so that F(pcrit)=infinity. This would suggest that the asymptotic behavior of the coefficients of the Maclaurin series for F(p) grows exponentially like (1/pcrit)^N more precisely their logs ought to be asymptotic. If this asymptotic could be proven with error bounds, then the problem of proving your conjecture that all coeffs are positive, would be reduced to a finite number of cases, which you could then handle by computer enumeration. There is a book Harry Kesten: Percolation, I have not looked at for 10 years, but might be a good start.