The first function (#1.) simply calculates the X^{2} statistic from the frequencies in variable y, and assumes y[1]=a, y[2]=b, y[3]=c, y[4]=d - where a,b,c,d, are the four cell frequencies in a 2×2 table. Note, the value of X^{2} is not available (NA) if any of the margin totals (a+b, or c+d, or a+c, or b+d) equal zero.
This function is needed for functions #3 and #4 to work.
The second function (#2.) is a random independent binomial function, and is needed in order to estimate the distribution of the X^{2} statistic assuming you have 2 samples, each of which represent binomially distributed populations.
The third function (#3.) takes R (=5000) random samples from a multinomial distribution, whose expected frequencies are calculated from the margin totals of your sample - assuming the null hypothesis to be true. The X^{2} statistic of each of these (bootstrap) samples is calculated, and arranged in ascending order. By finding what proportion of this bootstrap distribution equals or exceeds the observed value of X^{2} a (conventional or mid-P) 1-tailed P-value is returned.
The fourth function (#4.) does exactly the same thing, but uses the independent binomial function, and assumes that a+b=n1 and c+d=n2.