Correlation functions

CorrelationProfile

Structure that wraps the result of the intermittent_correlation_profile function: a set of intermittent correlation functions, one for each shell (bin) of distances to the binding site.

Fields

• r::Vector{T}: the center of each distance bin.
• correlations::Vector{OffsetVector{Float64,Vector{Float64}}}: for each bin, the intermittent correlation function, an OffsetArray with indices 0:maxdelta. Bins for which no solvent molecule was ever found contain NaN values.
• counts::Vector{Int}: for each bin, the number of (frame, solvent molecule) observations found in that bin. Useful to discard poorly sampled bins.
• delta_r::T: the width of each bin.
• step_r::T: the displacement between the lower edges of consecutive bins.
• dmax::T: the maximum distance considered, inherited from the input data.

intermittent_correlation(
    hbo::HydrogenBondOccupancy;
    maxdelta::Integer = length(hbo.list) ÷ 10,
    show_progress::Bool = true,
)

Calculate the intermittent correlation function of the hydrogen bonds found by the hydrogen_bond_occupancy function. That is, computes the probability of finding a given hydrogen bond (the same donnor, polar hydrogen, and acceptor atoms) present at frame i + delta, given that it was present at frame i, for each hydrogen bond independently. The hydrogen bond does not need to remain present in every frame in between: it may break and reform within the interval (that is what makes this correlation function "intermittent", as opposed to "continuous").

Returns an OffsetArray with indices 0:maxdelta, where the value at position 0 is 1.0, corresponding to the normalized count of events.

Arguments

• hbo::HydrogenBondOccupancy: The result of the hydrogen_bond_occupancy function.
• maxdelta::Integer: The maximum delta-step to be considered. Defaults to length(hbo.list) ÷ 10.
• show_progress::Bool: Show progress bar. Defaults to true.

Example

julia> using MolSimToolkit, MolSimToolkit.Testing

julia> sim = Simulation(Testing.namd_pdb, Testing.namd_traj);

julia> hbo = hydrogen_bond_occupancy(sim, "protein", show_progress=false)["protein => protein"];

julia> c = intermittent_correlation(hbo; maxdelta=4, show_progress=false);

julia> c[0]
1.0

intermittent_correlation(
    occupancy::Occupancy;
    maxdelta::Integer = length(occupancy.list) ÷ 10,
    show_progress::Bool = true,
)

Calculate the intermittent correlation function of the occupancy of a binding site, as computed by the occupancy function. That is, computes the probability of finding a solvent molecule at the site at frame i + delta, given that it was found at the site at frame i, for each solvent molecule independently.

Returns an OffsetArray with indices 0:maxdelta, where the value at position 0 is 1.0, corresponding to the normalized count of events.

Arguments

• occupancy::Occupancy: The result of the occupancy function.
• maxdelta::Integer: The maximum delta-step to be considered. Defaults to length(occupancy.list) ÷ 10.
• show_progress::Bool: Show progress bar. Defaults to true.

Example

julia> using MolSimToolkit, PDBTools, MolSimToolkit.Testing

julia> sim = Simulation(Testing.namd2_pdb, Testing.namd2_traj);

julia> protein = select(get_atoms(sim), "protein");

julia> tmao = select(get_atoms(sim), "resname TMAO");

julia> occ = occupancy(sim, protein, tmao; solvent_natomspermol=14, cutoff=3.0, show_progress=false);

julia> c = intermittent_correlation(occ; maxdelta=4, show_progress=false);

julia> c
5-element OffsetArray(::Vector{Float64}, 0:4) with eltype Float64 with indices 0:4:
 1.0
 0.39
 0.23469387755102042
 0.12903225806451613
 0.033707865168539325

intermittent_correlation(
    data::AbstractVector; 
    maxdelta = length(data) ÷ 10, 
    types::Function = x -> true,
    show_progress::Bool = true,
)

Calculate the intermittent correlation function of a time series. That is, computes the probability of finding a value of the same type at a step i + delta in the time series, given that it was present in step i.

Returns an OffsetArray with indices 0:maxdelta, where the value at position 0 is 1.0, corresponding to the normalized count of events.

Arguments

• data::AbstractVector: The time series to be analyzed.
• maxdelta::Integer: The maximum delta-step to be considered. Defaults to length(data) ÷ 10.
• types (optional): A function that returns true for the types of data that should be considered. Defaults to all data, i. e. x -> true. For example, to ignore 0 values, use types = x -> x != 0.
• show_progress::Bool: Show progress bar. Defaults to true.

Examples

Here we produce a time-series of 10,000 elements, as a sequence of 1's and 0's ([1, 0, 1, 0, ...]), and calculate the intermittent correlation function. The probability of finding the same number (0 or 1) after odd steps is 0, and the probability of finding the same number after even steps is 1.

julia> using MolSimToolkit

julia> data = [ mod(i,2) for i in 1:10^4 ];

julia> intermittent_correlation(data; maxdelta=4, show_progress=false)
5-element OffsetArray(::Vector{Float64}, 0:4) with eltype Float64 with indices 0:4:
 1.0
 0.0
 1.0
 0.0
 1.0

julia> intermittent_correlation(data; maxdelta=4, types = x -> x != 0, show_progress=false)
5-element OffsetArray(::Vector{Float64}, 0:4) with eltype Float64 with indices 0:4:
 1.0
 0.0
 1.0
 0.0
 1.0

In the second run, we have ignored the 0 values, and the result is the same, because here the correlations of the 1 values are the same as the correlations of the 0 values.

intermittent_correlation_profile(
    occ::Occupancy;
    delta_r::Real,
    step_r::Real = delta_r / 10,
    maxdelta::Integer = length(occ.list) ÷ 10,
    show_progress::Bool = true,
)

Computes the intermittent correlation function of the site occupancy resolved by the distance between the site and the solvent molecules, as computed by the occupancy function.

The distances from 0 to occ.dmax are split into bins of constant width delta_r, with the lower edge of consecutive bins displaced by step_r. Since step_r is typically smaller than delta_r, the bins overlap, and the profile is a smooth (quasi-continuous) function of the distance.

For each bin, the correlation function at delta is the probability that a solvent molecule found in that bin at frame i is found in the same bin at frame i + delta. The molecule does not need to remain in the bin in the frames in between: it may leave and come back (that is what makes this correlation function "intermittent"). Thus, each bin characterizes how long a solvent molecule remains at that distance from the site, and the decay of the correlation function can be converted into a residence time with the residence_time function.

Arguments

• occ::Occupancy: The result of the occupancy function. It must carry distance information, that is, it must have been produced by occupancy and not built from a list of molecule indices alone.

Keyword arguments

• delta_r::Real: The width of each distance bin. Required.
• step_r::Real: The displacement between the lower edges of consecutive bins. Defaults to delta_r / 10.
• maxdelta::Integer: The maximum delta-step to be considered. Defaults to length(occ.list) ÷ 10.
• show_progress::Bool: Show progress bar. Defaults to true.

Returns

• A CorrelationProfile object.

Example

julia> using MolSimToolkit, PDBTools, MolSimToolkit.Testing

julia> sim = Simulation(Testing.namd2_pdb, Testing.namd2_traj);

julia> protein = select(get_atoms(sim), "protein");

julia> tmao = select(get_atoms(sim), "resname TMAO");

julia> occ = occupancy(sim, protein, tmao; solvent_natomspermol=14, cutoff=6.0, show_progress=false);

julia> p = intermittent_correlation_profile(occ; delta_r=1.0, step_r=0.5, maxdelta=4, show_progress=false)
-------------------------------------------------------------------
Intermittent correlation profile:
-------------------------------------------------------------------
Number of distance bins: 11
Bin width (delta_r): 1.0
Bin step (step_r): 0.5
Maximum distance (dmax): 6.0
Range of bin centers: 0.5 - 5.5
Maximum delta (frames): 4
Observations per bin: minimum = 0, maximum = 89
-------------------------------------------------------------------

julia> p.r[end] # center of the outermost distance bin
5.5

julia> p.correlations[end][0:2] # correlation of the outermost bin
3-element Vector{Float64}:
 1.0
 0.05357142857142857
 0.019230769230769232

residence_time(c::AbstractVector; threshold::Real = 0.5, dt::Real = 1)
residence_time(p::CorrelationProfile; threshold::Real = 0.5, dt::Real = 1)

Returns the time at which an intermittent correlation function falls below threshold. For threshold = 0.5 (the default), this is the half-life of the correlation, that is, the time after which half of the solvent molecules initially present are no longer there.

When applied to a CorrelationProfile, returns one residence time for each distance bin, that is, the residence time as a function of the distance to the site.

The time is obtained by linear interpolation between the two consecutive delta-steps that bracket the threshold, and is thus not necessarily an integer number of frames. NaN is returned if the correlation function does not fall below threshold within the maxdelta of the input data (or if the bin contains no data at all).

Keyword arguments

• threshold::Real: The correlation value defining the residence time. Defaults to 0.5.
• dt::Real: The time interval between consecutive frames of the input data, used to convert the result from frames to time units. Defaults to 1, that is, the result is given in frames.

Example

julia> using MolSimToolkit, PDBTools, MolSimToolkit.Testing

julia> sim = Simulation(Testing.namd2_pdb, Testing.namd2_traj);

julia> protein = select(get_atoms(sim), "protein");

julia> tmao = select(get_atoms(sim), "resname TMAO");

julia> occ = occupancy(sim, protein, tmao; solvent_natomspermol=14, cutoff=6.0, show_progress=false);

julia> ic = intermittent_correlation(occ; maxdelta=4, show_progress=false);

julia> residence_time(ic) # in frames
0.9035714285714286

julia> p = intermittent_correlation_profile(occ; delta_r=1.0, step_r=0.5, maxdelta=4, show_progress=false);

julia> residence_time(p)[end] # outermost distance bin
0.5283018867924528