OUwie (version 3.0.3)

hOUwie.walk: Sample points from along a ridge for a hOUwie model

Description

This "dents" the likelihood surface by reflecting points better than a threshold back across the threshold (think of taking a hollow plastic model of a mountain and punching the top so it's a volcano). It then uses essentially a Metropolis-Hastings walk to wander around the new rim. It adjusts the proposal width so that it samples points around the desired likelihood. This is better than using the curvature at the maximum likelihood estimate since it can actually sample points in case the assumptions of the curvature method do not hold. It is better than varying one parameter at a time while holding others constant because that could miss ridges: if I am fitting 5=x+y, and get a point estimate of (3,2), the reality is that there are an infinite range of values of x and y that will sum to 5, but if I hold x constant it looks like y is estimated very precisely. Of course, one could just fully embrace the Metropolis-Hastings lifestyle and use a full Bayesian approach.

While running, it will display the current range of likelihoods in the desired range (by default, the best negative log likelihood + 2 negative log likelihood units) and the parameter values falling in that range. If things are working well, the range of values will stabilize during a search.

Usage

hOUwie.walk(houwie_obj, 
            delta                 = 2, 
            nsteps                = 1000, 
            print_freq            = 50, 
            lower_bound           = 0, 
            upper_bound           = Inf, 
            adjust_width_interval = 100, 
            badval                = 1e9, 
            sd_vector             = NULL, 
            debug                 = FALSE, 
            restart_after         = 50)

Value

A dentist object containing results, the data.frame of negative log likelihoods and the parameters associated with them; acceptances, the vector of whether a proposed move was accepted each step; best_neglnL, the best value passed into the analysis; delta, the desired offset; all_ranges, a summary of the results.

Arguments

houwie_obj

A list of class "houwie".

delta

How far from the optimal negative log likelihood to focus samples.

nsteps

How many steps to take in the analysis.

print_freq

Output progress every print_freq steps.

lower_bound

Minimum parameter values to try. One for all or a vector of the length of par.

upper_bound

Maximum parameter values to try. One for all or a vector of the length of par.

adjust_width_interval

When to try automatically adjusting proposal widths.

badval

Bad negative log likelihood to return if a non-finite likelihood is returned.

sd_vector

Vector of the standard deviations to use for proposals. Generated automatically if NULL.

debug

If TRUE, prints out much more information during a run.

restart_after

Sometimes the search can get stuck outside the good region but still accept moves. After this many steps without being inside the good region, restart from one of the past good points.

Details

The algorithm tunes: if it is moving too far away from the desired likelihoods, it will decrease the proposal width; if it staying in areas better than the desired likelihood, it will increase the proposal width. It will also expand the proposal width for parameters where the extreme values still appear good enough to try to find out the full range for these values.

In general, the idea of this is not to give you a pleasingly narrow range of possible values -- it is to try to find the actual uncertainty, including finding any ridges that would not be seen in univariate space.

Examples

# ## Finally, it is important to get estimates of your parameter uncertainty. 
# ## For that we have imported the functions of Brian O'Meara's dentist package
# data(tworegime)
# 
# 
# ##Third step is to fit our model (here we are using fixed params from a previous ML search):
# p <- c(0.01664314, 0.39631218, 0.18684476, 2.25121568, 0.82495093) # MLE
# pp_cd <- hOUwie(tree, trait, rate.cat = 1, discrete_model = "ER", 
#                 continuous_model = "OUM", nSim = 25, p = p)
# 
# ## We can compare this to the result from OUwie on a fixed map to see a difference 
# dent_res <- hOUwie.walk(pp_cd, nsteps = 25, debug =TRUE)
# plot(dent_res)