S^2 parameter in gmx chi

User discussions
1

GROMACS version:5.1.2
GROMACS modification: No
Here post your question

Dear Gromacs User,

I have come across the gmx chi and was wondering which equation does gromacs use to calculate the S^2 parameter (in the order.xvg). Furthermore, by follow the paper recommended by the software, “Protein chemical shift analysis: a practical guide”, the delta Ca potential is it calculate based on their formula and from the formula, I am not too sure how random coil chemical shifts fit into that? I would appreciate if someone can point to where those equations can be found and where it came from. Thank you very much.

Best regards,

Ben

2

Hi Ben,

I’m quoting the part of the source code where the actual evaluation is happening in calc_order in gmx_order.cpp below and post what I would translate this to in terms of math with a direction \vec{v}

\vec{S}_z = \cfrac{\vec{C}_n -\vec{C}_{n+1}}{ \|\vec{C}_n -\vec{C}_{n+1}\|} ; \vec{S}_x = (\vec{C}_n -\vec{C}_{n-1} )\times(\vec{C}_n -\vec{C}_{n+1} ); \vec{S}_y = \vec{S}_x \times \vec{S}_z;

S^2 = \left\langle\frac32 \left((\vec{S}_x\cdot \vec{v})^2 +(\vec{S}_y\cdot \vec{v})^2+(\vec{S}_z\cdot \vec{v})^2 \right) -\frac12 \right\rangle_{\mathrm{atoms}}

            /* now get Sx. Sx is normal to the plane of Cn-1, Cn and Cn+1 so
               we can use the outer product of Cn-1->Cn and Cn+1->Cn, I hope */
            rvec_sub(x1[a[index[i + 1] + j]], x1[a[index[i] + j]], tmp1);
            rvec_sub(x1[a[index[i - 1] + j]], x1[a[index[i] + j]], tmp2);
            cprod(tmp1, tmp2, Sx);
            svmul(1.0 / norm(Sx), Sx, Sx);

            /* now we can get Sy from the outer product of Sx and Sz   */
            cprod(Sz, Sx, Sy);
            svmul(1.0 / norm(Sy), Sy, Sy);

            /* the square of cosine of the angle between dist and the axis.
               Using the innerproduct, but two of the three elements are zero
               Determine the sum of the orderparameter of all atoms in group
             */
            if (use_unitvector)
            {
                cossum[XX] = gmx::square(iprod(Sx, direction)); /* this is allowed, since Sa is normalized */
                cossum[YY] = gmx::square(iprod(Sy, direction));
                cossum[ZZ] = gmx::square(iprod(Sz, direction));
            }
            else
            {
                cossum[XX] = gmx::square(Sx[axis]); /* this is allowed, since Sa is normalized */
                cossum[YY] = gmx::square(Sy[axis]);
                cossum[ZZ] = gmx::square(Sz[axis]);
            }

            for (m = 0; m < DIM; m++)
            {
                frameorder[m] += 0.5 * (3.0 * cossum[m] - 1.0);
            }
       for (m = 0; m < DIM; m++)
        {
            (*order)[i][m] += (frameorder[m] / static_cast<real>(size));
        }
3

Hi cblau,

Thank you very much. Can you share the path to that file? I tried to look for it and can’t find it. Thank you very much.

Best,

Ben

4

Hi Ben,

I checked what gmx order does instead of gmx chi!
The essential part for gmx chi is here:

Sorry for the confusion!