- Contemporary messages sorted: [ by date ] [ by thread ] [ by subject ] [ by author ] [ by messages with attachments ]

From: Osman, Roman <roman.osman.mssm.edu>

Date: Thu, 1 Sep 2016 16:06:44 +0000

Hi Tom and anyone who can help,

After some experimentation with the cpptraj autocorr and our own code we think that we have reconciled the differences between the two calculations.

I would be grateful if anyone can corroborate our observation and comment whether this is indeed what the code in amber is doing.

We noticed that short tracks of h-bonds that go rapidly to 0, fit very well to those that we calculate with our code. But tracks that decay over a long period of time, i.e., longer than about 300 ns, the auc calculated with amber had a larger slope with essentially the same shape as our auc. This led to the supposition that amber may be adding zeros at the end of the track to ensure that the limit of the track is going to zero at infinite time.

Indeed, when we added zeros at the end of the tracks and used our program the auc were identical.

Of course, that is not a proof that this is what is done in cpptraj.

So we would like to get a possible answer of how autocorr in cpptraj is doing the calculation.

As I wrote before we use the nocovar method because otherwise the auc becomes negative at long times.

Thanks,

Rami

On Aug 30, 2016, at 6:07 PM, Tom Kurtzman <simpleliquid.gmail.com<mailto:simpleliquid.gmail.com>> wrote:

Hi Roman, in response to your comment "In my opinion this is not possible

since the track is made of 0 and 1, so the autocorrelation cannot be < 0."

Often correlation functions (including autocorrelation functions) are

normalized so that their initial values are 1 and absence of correlation is

zero so instead of calculating

<A(t)A(t')> they actually calculate This is the same as calculating

<[A(t) -<A>][A(t')-<A>]> = <A(t)A(t')> - <A>^2. This allows for negative

numbers and ensures that the limit as t --> \infty is zero. Often the

function is then divided by its initial value to give an initial value

(upper bound) value of unity.

For the no-covariance calculation, I assume that the average values are not

subtracted and hence the limit and infinite time approaches <A>^2 which is

positive.

It is not unusual for autocorrelation functions to have negative values nor

to have long time non-zero tails (the velocity autocorrelation function is

a famous example of both of these) however, I'm not sure whether the long

time tail you are seeing is

1) real

2) fluctuations about zero due to sampling limitations (noise)

3) an artifact of an incorrectly calculated quantity.

Without seeing the data, I would initially suspect it's just noise.

Tom

On Mon, Aug 29, 2016 at 9:05 AM, Osman, Roman <roman.osman.mssm.edu<mailto:roman.osman.mssm.edu>> wrote:

Hello Amberites,

I was wondering whether anybody has used an autocorrelation of hydrogen

bonds.

Here is a cpptraj input that I ran:

# Read in trajectory

trajin mcmrjp_free_nw.dcd

hbond HB out mcmrjp_free_hb.dat :17,252,67,206,159,272,27,32,

228,245,262,267,165

,183,34,138,239,290,39,143,207,204,162 series uuseries mcmrjp_free_hb.dat

avgout

mcmrjp_free_avg.dat printatomnum

run

runanalysis lifetime HB[solutehb] out mcmrjp_free_lifet.dat

autocorr HB[solutehb] out mcmrjp_free_auc.dat lagmax 14500

run

The surprising thing is that I get in the aucorrelation function negative

numbers towards the end of the track.

In my opinion this is not possible since the track is made of 0 and 1, so

the autocorrelation cannot be < 0.

In fact if I run with

autocorr HB[solutehb] out mcmrjp_free_auc.dat lagmax 14500 nocovar

The autocorrelations become positive.

Am I doing something wrong?

Thanks for your help.

PS: I can provide the trajectory (large) and the output from the

autocorrelation.

Roman Osman

Professor of Structural and Chemical Biology

Mount Sinai School of Medicine

New York, NY 10029

(212) 659-8627

roman.osman.mssm.edu<mailto:roman.osman.mssm.edu><mailto:roman.osman.mssm.edu>

_______________________________________________

AMBER mailing list

AMBER.ambermd.org<mailto:AMBER.ambermd.org>

https://urldefense.proofpoint.com/v2/url?u=http-3A__lists.ambermd.org_mailman_listinfo_amber&d=AwICAg&c=4R1YgkJNMyVWjMjneTwN5tJRn8m8VqTSNCjYLg1wNX4&r=80h8gfWnXRd-mXPSrw98JBRjhncHe1Sopk2U91C1sZg&m=-hoofXSvixmoA1iOwfXi41mSpf3DaEZs1kdxf6s0p-Q&s=YXoZ0b2W1iVWNm4t_-F4rjkMCoaEH-OAP15eidtQO2s&e=

Date: Thu, 1 Sep 2016 16:06:44 +0000

Hi Tom and anyone who can help,

After some experimentation with the cpptraj autocorr and our own code we think that we have reconciled the differences between the two calculations.

I would be grateful if anyone can corroborate our observation and comment whether this is indeed what the code in amber is doing.

We noticed that short tracks of h-bonds that go rapidly to 0, fit very well to those that we calculate with our code. But tracks that decay over a long period of time, i.e., longer than about 300 ns, the auc calculated with amber had a larger slope with essentially the same shape as our auc. This led to the supposition that amber may be adding zeros at the end of the track to ensure that the limit of the track is going to zero at infinite time.

Indeed, when we added zeros at the end of the tracks and used our program the auc were identical.

Of course, that is not a proof that this is what is done in cpptraj.

So we would like to get a possible answer of how autocorr in cpptraj is doing the calculation.

As I wrote before we use the nocovar method because otherwise the auc becomes negative at long times.

Thanks,

Rami

On Aug 30, 2016, at 6:07 PM, Tom Kurtzman <simpleliquid.gmail.com<mailto:simpleliquid.gmail.com>> wrote:

Hi Roman, in response to your comment "In my opinion this is not possible

since the track is made of 0 and 1, so the autocorrelation cannot be < 0."

Often correlation functions (including autocorrelation functions) are

normalized so that their initial values are 1 and absence of correlation is

zero so instead of calculating

<A(t)A(t')> they actually calculate This is the same as calculating

<[A(t) -<A>][A(t')-<A>]> = <A(t)A(t')> - <A>^2. This allows for negative

numbers and ensures that the limit as t --> \infty is zero. Often the

function is then divided by its initial value to give an initial value

(upper bound) value of unity.

For the no-covariance calculation, I assume that the average values are not

subtracted and hence the limit and infinite time approaches <A>^2 which is

positive.

It is not unusual for autocorrelation functions to have negative values nor

to have long time non-zero tails (the velocity autocorrelation function is

a famous example of both of these) however, I'm not sure whether the long

time tail you are seeing is

1) real

2) fluctuations about zero due to sampling limitations (noise)

3) an artifact of an incorrectly calculated quantity.

Without seeing the data, I would initially suspect it's just noise.

Tom

On Mon, Aug 29, 2016 at 9:05 AM, Osman, Roman <roman.osman.mssm.edu<mailto:roman.osman.mssm.edu>> wrote:

Hello Amberites,

I was wondering whether anybody has used an autocorrelation of hydrogen

bonds.

Here is a cpptraj input that I ran:

# Read in trajectory

trajin mcmrjp_free_nw.dcd

hbond HB out mcmrjp_free_hb.dat :17,252,67,206,159,272,27,32,

228,245,262,267,165

,183,34,138,239,290,39,143,207,204,162 series uuseries mcmrjp_free_hb.dat

avgout

mcmrjp_free_avg.dat printatomnum

run

runanalysis lifetime HB[solutehb] out mcmrjp_free_lifet.dat

autocorr HB[solutehb] out mcmrjp_free_auc.dat lagmax 14500

run

The surprising thing is that I get in the aucorrelation function negative

numbers towards the end of the track.

In my opinion this is not possible since the track is made of 0 and 1, so

the autocorrelation cannot be < 0.

In fact if I run with

autocorr HB[solutehb] out mcmrjp_free_auc.dat lagmax 14500 nocovar

The autocorrelations become positive.

Am I doing something wrong?

Thanks for your help.

PS: I can provide the trajectory (large) and the output from the

autocorrelation.

Roman Osman

Professor of Structural and Chemical Biology

Mount Sinai School of Medicine

New York, NY 10029

(212) 659-8627

roman.osman.mssm.edu<mailto:roman.osman.mssm.edu><mailto:roman.osman.mssm.edu>

_______________________________________________

AMBER mailing list

AMBER.ambermd.org<mailto:AMBER.ambermd.org>

https://urldefense.proofpoint.com/v2/url?u=http-3A__lists.ambermd.org_mailman_listinfo_amber&d=AwICAg&c=4R1YgkJNMyVWjMjneTwN5tJRn8m8VqTSNCjYLg1wNX4&r=80h8gfWnXRd-mXPSrw98JBRjhncHe1Sopk2U91C1sZg&m=-hoofXSvixmoA1iOwfXi41mSpf3DaEZs1kdxf6s0p-Q&s=YXoZ0b2W1iVWNm4t_-F4rjkMCoaEH-OAP15eidtQO2s&e=

-- ************************************************ Tom Kurtzman, Ph.D. Assistant Professor Department of Chemistry Lehman College, CUNY 250 Bedford Park Blvd. West Bronx, New York 10468 718-960-8832 https://urldefense.proofpoint.com/v2/url?u=http-3A__www.lehman.edu_faculty_tkurtzman_&d=AwICAg&c=4R1YgkJNMyVWjMjneTwN5tJRn8m8VqTSNCjYLg1wNX4&r=80h8gfWnXRd-mXPSrw98JBRjhncHe1Sopk2U91C1sZg&m=-hoofXSvixmoA1iOwfXi41mSpf3DaEZs1kdxf6s0p-Q&s=o3uJDn8But0uJdUcDMpoKaCqj2NRM6dxD3CbDdwWrHM&e= <https://urldefense.proofpoint.com/v2/url?u=http-3A__www.lehman.edu_faculty_tkurtzman_index.html&d=AwICAg&c=4R1YgkJNMyVWjMjneTwN5tJRn8m8VqTSNCjYLg1wNX4&r=80h8gfWnXRd-mXPSrw98JBRjhncHe1Sopk2U91C1sZg&m=-hoofXSvixmoA1iOwfXi41mSpf3DaEZs1kdxf6s0p-Q&s=G9qKDa1ESYL7tKFnPsa8zDrGteWJqxHWTUzqap0HM78&e= > ************************************************ _______________________________________________ AMBER mailing list AMBER.ambermd.org<mailto:AMBER.ambermd.org> https://urldefense.proofpoint.com/v2/url?u=http-3A__lists.ambermd.org_mailman_listinfo_amber&d=AwICAg&c=4R1YgkJNMyVWjMjneTwN5tJRn8m8VqTSNCjYLg1wNX4&r=80h8gfWnXRd-mXPSrw98JBRjhncHe1Sopk2U91C1sZg&m=-hoofXSvixmoA1iOwfXi41mSpf3DaEZs1kdxf6s0p-Q&s=YXoZ0b2W1iVWNm4t_-F4rjkMCoaEH-OAP15eidtQO2s&e= Roman Osman Professor of Structural and Chemical Biology Mount Sinai School of Medicine New York, NY 10029 (212) 659-8627 roman.osman.mssm.edu<mailto:roman.osman.mssm.edu> _______________________________________________ AMBER mailing list AMBER.ambermd.org http://lists.ambermd.org/mailman/listinfo/amberReceived on Thu Sep 01 2016 - 09:30:02 PDT

Custom Search