WeiYa's Work Yard

A dog, who fell into the ocean of statistics, tries to write down his ideas and notes to save himself.

Growing A Polymer

Posted on July 17, 2017 0 Comments

Model: Self-avoid walk

For the sake of a concise description of the model, we just consider two-dimensional lattice space.

In such a model, the realization of a chain polymer of length $N$ is fully characterized by the positions of all of its molecules, $\mathbf x = (x_1, x_2,\ldots, x_N)$, where $x_i = (a, b)$ is a point in the 2-D lattice space. The distance between $x_i$ and $x_{i+1}$ has to be exactly 1, and $x_{i+1}\neq x_k, k<i$

A possible Boltzmann distribution for the chain polymer modeled by an SAW is the uniform distribution. More formally, the target distribution of interest in this case is where the space of $\mathbf x$ is the set of all SAWs of length $N$ and $Z_N$ is the normalizing constant, which is just the total number of different SAWs with $N$ atoms.

Let $x_1 = (0, 0), x_2=(1, 0)$. Draw the position of $x_{t+1}$ conditional on the current configuration of $(x_1, \ldots, x_t)$, according to the probability distribution

where $(i’, j’)$ is one of the unoccupied neighbors of $x_t$ and $n_t$ is the total number of such unoccupied neighbors.

Unfortunately, the SAWs produced by this way is not uniformly distributed. Considering the target distribution is $\pi(\mathbf x)\propto 1$ and the sampling distribution of $\mathbf x$ is $(n_1\times \cdots\times n_{N-1})^{-1}$, so assign a weigh

to the SAW.

Implement this simple model by using the following R code

Sequential importance sampling

Apply the sequential importance sampling, there are two alternatives for auxiliary sequence of distribution.

Firstly, consider the sequence of auxiliary distributions $\pi_t(\mathbf x_t), t=1,\ldots, N-1$ is just the sequence of auxiliary distributions of SAWs with $t$ nodes, where $Z_t$ is the total number of such SAWs. Under $\pi_t$, the marginal distribution of the partial sample $\mathbf x_{t-1}=(x_1,\ldots, x_{t-1})$ is

Then, $\pi_t(x_t\mid \mathbf x_{t-1})=1/n_{t-1}$, so we can choose the sampling distribution $g_t(x_t\mid \mathbf x_{t-1})$ as $1/n_{t-1}$. If $n_{t-1} = 0$, $g(x_t\mid \mathbf x_{t-1}) = 0.$

In addition to this auxiliary distributions, we can use the marginal distribution of $\mathbf x_t$ under the uniform distribution of $(\mathbf x_t, x_{t+1})$

then we can let $g_t$ be $\pi_t’(x_t\mid \mathbf x_{t-1})$. This procedure is equivalent to a two-step-look-ahead approach.

We can compute certain descriptive statistics regarding such SAWs, say the mean squared extension of the chain $E\Vert x_N-x_1\Vert^2$ or the total number of possible configuration of SAWs with $N$ nodes $Z_N$.

In order to estimate $Z_N$, consider the importance weight. We can compute the importance weight recursively as $w_{t+1}=w_tn_t$. Then we obtain

(??: Is $Z_t/Z_{t+1}$ const? And is it ok to omit this term?)

The target distribution is $\pi(\mathbf x) = \frac{1}{Z_N}$, the weight satisfies the relationship

Hence, $E_p(w_N) = Z_N$, which gives us a means to estimate the partition function (normalizing constant):

where $\bar w_N$ is the sample average of the weights obtained by repeating the biased sampling procedure multiple $m$ times.

References

Liu, Jun S. Monte Carlo strategies in scientific computing. Springer Science & Business Media, 2008.


Published in categories Monte Carlo