Draft
Conversation
First compiling version of major changes; still optimisation and further structural changes yet to do.
Originally, SIFT2 had two different modes of regularisation operating in parallel. These were termed "Tikhonov" ("tik" for short) and "asymmetric Total Variation" ("aTV" for short). Firstly, the name of the latter was erroneous. Secondly, this did not properly cover the space of possibilities, and was consequently restrictive of future prospects.
These changes offer two "bases" of regularisation, with three "functions" available for each. But unlike the original implementation, only one form of regularisation can be used in any single invocation.
- "Basis": This is where regularisation is applied either to each streamline entirely independently of all others, or based on the difference between the contribution of a streamline and the other streamlines traversing each fixel it traverses.
- "Function": This is the shape of the function that controls how variance in streamline contributions are penalised. The "gamma" function proposed in the original manuscript, which penalises by coefficient below the mean and weight above the mean, is one of the three options.
The two regularisation forms offered in the original implementation are therefore two of six possible options in this new design.
Hopefully fixes compilation issues as reported in #2674.
Note that this also includes modifications to the so-called "gaussian" branch of the mapper, which may improve calculations for streamlines where vertices are not uniformly spaced and-or if used in conjunction with -precise.
Need ability to distinguish between the header that defines the output image, and the header that defines the voxel grid to which the streamlines are mapped, since for fixels these are drastically different. While existing variable "template_header" theoretically could have been re-used for this purpose, I chose to leave that fixed as what was read from the input to the -template option if used and invalid if not used.
The SIFT2 algorithm is augmented to support attribution of both absolute fixel-wise fibre densities, and fixel-wise fibre density differences, to streamlines. - Addition of options to initialise or set absolute streamline weights prior to running in differential mode. - Better conformance of referring to absolute weights as "coefficient" (what appears in the exponent) and "factor" (what scales the contribution of the streamline in the model). - Addition of a novel absolute weight optimiser based on Barzilai-Borwein gradient descent; this yields unstable behaviour, and has not yet been debugged. - New class WeightingCoeffAndFactor is used to reduce redundant std::log() and std::exp() calls by immediately quantifying both weighting forms when a streamline is to be processed. This is also used to provide better disambiguation between functions intended to be used for absolute vs. differential modes.
Note also slight behavioural change to -in_factor / -in_coeffs: Limits on minimum / maximum coefficients / factors in absolute mode will not be imposed given that no absolute-mode optimisation takes place if either of these two options is specified.
New members within the TckFactor class provide a binary mask across streamlines (one for absolute mode, one for differentiail mode), and these control whether or not a streamline can contribute to the optimisation in the next iteration.
The documentation indicates that if a streamline obtains a weighting coefficient that is equal to or lesser than the minimum permissible, then it should be removed from the optimisation altogether; that is, given a weighting coefficient of -inf and a weighting factor of 0. However looking at the code it seems that this was not being applied.
- Disable fixel-based regularisation in differential mode until new expressions are derived. - Update data-driven cost function and derivatives for differential mode under original formulation. - Introduce a tailored regularisation function "asymptotic" for differential mode.
Primary output from differential mode is now the product of the absolute-mode weighting factor and the delta coefficient multiplier, therefore representing a FBC change.
This commit is an amalgamation of multiple changes that were initially applied to the experimental "asymmetric" SIFT2 differential, back-propagated to the original SIFT2 differential formulation. This includes: - Fixing command operation when streamline differential weights are specified as fixed vs. initialised. - Improved code and documentation regarding the definitions of "delta coefficient" versus "delta weight" (note that the definitions of these terms differ between this commit and its origin as developments for the "SIFT2 asymmetric differential" mode; similar clarification has been made here but based on the definions of the original formulation). This includes new class "DifferentialWCF", which forces code to be explicit about which variable is to be used. - Shuffling of regularisation code, so that functions relevant to absolute mode are always presented in the same order (coeffieitn, factor, gamma). - Handle NaN input FD values (which presumably arises if utilising the FDC metric per fixel in a template but FC attempts to sample from outside the warp image FoV). - Set termination of optimiser based on a fraction of the naive cost function for absolute mode, rather than its initial value, which may be influenced by initial values or be reduced in magnitude due to two sessions in differential mode being closer to one another than to the naive reconstruction. - Restore calculation of fixel-wise mean delta coefficients in differential mode. - Fix manifestation of NaN values where a streamline does not contribute to any fixel with a non-zero model weight. - Ensure that manually-specified numbers of iterations apply to both absolute and differential modes. - New option -in_mu, which allows setting a fixed value of mu from some other source that is irrespective of the FD / TD data for the active reconstruction; this may be of particular utility when providing initial weighting coefficients / factors from another reconstruction and wishing for them to retain their original scaling.
Derived and implemented derivative calculations for new differential-mode-specific regularisation functions when operating on a per-fixel basis.
- Constrain domain of delta coefficient depending on functional form of regulariser. - Increase requisite precision of solution per line search.
…sed regularisation
- Increase precision of line search. - Preclude streamlines with a delta coefficient of -1.0 or 1.0 (which should only be possible via user initialisation) from sending the total regularisation cost to infinity. - Fix algebraic error in dual inverse barrier function for fixel-basis regularisation. - Yield an error in the scenario where the user initialises delta coefficients to +- 1.0, but then fails to explicitly exclude those same streamlines from the optimisation, if the dual inverse barrier function is active for regularisation, as its derivatives do not exist there.
….0; previously checked if absolute deltacoeff is >= than -1, which always was true preventing execution of the command. now checks if absolute deltacoeff is >= than 1
- Prevent application of differential-mode coefficient restrictions to absolute-mode data. - Improve calculation of number of streamlines participating in optimisation in first iteration.
The -blur_streamlines option performs something akin to trilinear interpolation when mapping streamlines to a voxel grid. Internally it uses a very similar mechaism to that invoked by the -precise option. It is currently available for tckmap, tcksift, tcksift2, fixelconnectivity.
This was sometimes resulting in dereferencing a null pointer. Hopefully adding it within usage(), rather than concatenating it to a const App::OptionGroup, will fix. Follows 8225a02.
Error out appropriately if streamline-group-wise regularisation is requested, but no streamline grouping information is provided.
Implements more explicit checking of construction values for class WeightingCoeffAndFactor.
- Modify detection of violation of constraints. In absolute mode such detections are made in an outer loop, whereas for differential mode these constraints are more intrinsic and do not want to expose to floating-point precision issues. - When using the dual inverse barrier regularisation function, use Halley's method rather than Newton-Raphson updates for finding the minimum, as that function is strongly non-quadratic. - Where a polynomial update results in cost function increase, perform a bisection search to find a smaller step in that direction that results in a decrease. This prevents the regularisation term of the cost function from blowing up due to landing near an inverse barrier. - Subtle change to dual inverse barrier regularisation functions: a value of exactly +- 1.0 will yield infinite cost and NaN derivatives, whereas an absolute value greater than 1 will yield NaNs across the board. - Fix computation of fixel-wise mean delta coefficients. Inconsistency in whether absolute mode factors should influence this result could cause fixel means that were outside of the range [-1.0, 1.0].
- For fixel-wise and group-wise regularisation bases, explicitly compute the denominator of the mean calculation rather than relying on pre-computation of the unmodulated absolute track density, as some streamlines may need to be excluded from that calculation. - For reporting of total regularisation cost, omit from the calculation those streamlines that are explicitly deemed to be absent from the model; whether by a weighting factor of 0 / weighting coefficient of -inf, or with a delta coefficient of -1 and exclusion from the differential mask.
If all streamlines ascribed to a fixel / streamline group attain a delta coefficient of +-1.0, and that becomes the reference value for that fixel / streamline group, then the cost function for regularising by delta coefficient with respect to that reference experiences a discontinuity. This change should yield sensible values when that occurs.
github-project-automation
bot
moved this to Backlog
in Longitudinal quantitative streamline tractography
|
Reference |
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Software development work to support an imminent publication from @ppruc.
Draft Pull Request: there is still considerable work to do to resolve these changes against the current
devtip, and other software commands simply will not compile with the current changeset as they have not yet been updated to make use of the new back-end.Summary
This augments the
tcksift2command to have two modes of operation:The historical SIFT2 method. Streamline weights are determined that converge the tractogram toward fixel-wise measures of absolute fibre density.
The quantitative data utilised per fixel is proportional not to the total fibre density estimated for that fixel, but to the fibre density change associated with some other experimental variable (eg. time for longitudinal analysis).
Note that the software interface has changed slightly, including for the pre-existing Absolute mode. The input is no longer an FOD image (to which the FOD segmentation algorithm was performed internally), but a fixel data file. One of the benefits of performing the FOD segmentation internally was the preservation of dixel information, which enabled assignment of streamlines to fixels taking into account FOD lobe shape, rather than merely the dot product with the fixel mean direction. This functionality is however preserved by saving that information (#2625) and loading it in conjunction with the fixel dataset (#2657).
Notes for any early adopters
You will need to compile from source (perhaps as a container to avoid MRtrix version clashing).
tcksift2command.Other commands in the software are currently broken by these changes.
In its current state this is achieved by running "
./build bin/tcksift2"; but these instructions will change as this code branch is brought up to date.Still to do
(not guaranteed to be complete at time of posting)
devtipThis will include filesystem restructuring,
clang-format, anything thatcheck_syntaxcomplains about.constexprs.@ppruc: I know we said change to streamline basis and leave at 0.1, but it occurred to me that I had been concerned about over-regularisation at the current defaults, so there's some concern that going to streamline basis, which can be at greater danger of precluding correction of the largest discrepancies, may go too hard. Might need to generate a few to convince myself it's not a problem.
Differential mode expects as input half-differences, for the delta coefficients to be bounded as (-1.0, +1.0). However I think it's just trusting the user to do so. One option maybe worth adding is that if the differential FD is more than half the absolute FD in any fixel, that's likely due to a user providing full-differences rather than half-differences.