Linear mixed model in R; modelling fixed effects with multiple levels and interactions. Help! I am new to R and to mixed linear modelling. I have a dataset with variables from a cross-sectional study looking at fractional anisotrophy (a property of the brains white matter) in 6 different white matter fibre tracts in the brain. For each fibre tract there are 2 measures (one from each hemisphere). There are 66 participants divided into two groups. We want to control the group comparisons for differences in age, the average FA volume across the brain (wholebrain FA) and tract volume.
I assume that Group (Patients/Controls), Tract (CB/SLF1/SLF2/SLF3/UF/OFST) and Hemisphere (Left/Right) are fixed effects and that Subject (n=66) is a random effect. I also assume that Age, Wholebrain FA and Tract volume should be modelled as fixed effects. For Age and Wholebrain FA there is one value for each participant, but for Volume there is one value for each observation. The attached picture presents the table in the long format. For each subject there are 12 observations.

We hypothesized that there would be group difference in each of the six fibre tracts. We had no a priori assumptions about the hemispheres but would like to explore this post-hoc. We would also like to explore associations between Age and FA in different tracts.
My suggested model look like this
mixed.lmer <- lmer(FA ~ Age + Wholebrain_FA + Volume + GroupTractHemisphere + (1|Subjects), data = DTI)

Question 1: Given that Tract and Hemisphere are assumed to be fixed variables but also are within-subject variables, are they correctly modelled? I am having a hard time understanding how the model "understand" that these variables have multiple levels from the way it is written above.
Question 2: The Volume variable is a within-subject variable whereas Age is a between subject variable. Should they then not be modelled differently?
Question 3: Whether or not to include a three-way interaction is a major debate in my research group. Some saying that for practical purposes its impossible to really make sense of it. Other say it can guide the decision of whether or not to test differences between groups for each hemisphere in each tract. Including a three way interaction to the model is likely to change the results significantly so it seems pretty important to get it right the first time. Any thought on this? Is it being a criminal to include it?
 A: 
Question 1: Given that Tract and Hemisphere are assumed to be fixed variables but also are within-subject variables, are they correctly modelled? I am having a hard time understanding how the model "understand" that these variables have multiple levels from the way its written above.

In most software, such as lme4 or GLMMadaptive it is not necessary to specify at which level a variable varies because, contrary to your understanding, the software really does "know". The level at which a variable varies is a property of the data and it is easy  to demonstrate with cross-tabulations. 
You may also want to allow a within-subject fixed effect to vary randomly across subjects in which case you can also specify it as a random slope. For example: 
lmer(FA ~ Age + Wholebrain_FA + Volume + Hemisphere + (Hemisphere | Subjects)

will estimate a fixed effect for Hemisphere and also allow it to vary by subject. The software will estimate a variance for the "random slope" of Hemisphere. 
The difference between the model without random slopes and with random slopes is that in the former, the "within-subject" variable is estimated to have a fixed effect which is the same for all subjects, whereas fitting random slopes allows each subject to have it's own effect of that variable (a global fixed effect and a random offset) 

Question 2: The Volume variable is a within-subject variable whereas Age is a between subject variable. Should they then not be modelled differently?

Fixed effects are estimated in the same way regardless of whether they vary within levels of a grouping variable (Subject in your case). This means that the entries in the model matrix of fixed effects will be quite different for within vs. between variables, but this is not something you need to worry about. These kinds of concerns often arise when people come from a traditional ANOVA background. 

Question 3: Whether or not to include a three-way interaction is a major debate in my research group. Some saying that for practical purposes its impossible to really make sense of it. Other say it can guide the decision of whether or not to test differences between groups for each hemisphere in each tract. Including a three way interaction to the model is likely to change the results significantly so it seems pretty important to get it right the first time. Any thought on this? Is it being a criminal to include it?

In general there is no problem in interpreting statistical interactions. They have a fairly simple interpretation. This question is too broad to answer. I would suggest posting a new question about this, and including as much detail as possible.
A: You are currently only including a random intercept, that is the base rate of $FA$ can differ between subjects. While that is a good first step for a multilevel model, you should usually also include random intercepts, in other words, the effect sizes should be allowed to differ between subjects. This closely relates to your question about within-subject variables.


*

*Your model understands that the categorical variables have multiple levels, because these variables can't be interpreted numerically directly. What number would correspond to "left" or "right"? However, multilevel models (basically an advanced form of regression analysis) always require numerical predictors. So whenever R encounters a categorical variable, it automatically recodes it into numerical values using dummy coding. Dummy coding is used because in most cases it is not too wrong. If you need more control over the coding (e.g. choosing the reference class of the dummy coding, centering, ANOVA coding, etc), there are additional methods that you may use.

*Categorial variables (with multiple levels) and within-subject variables are two quite different concepts you seem to confuse somewhat. If you had included handedness or gender, those would have been between-subject but categorical (with different levels instead of numerically). It's good to keep those two concepts separate.

*The (1|Subjects) term in your model indicates that the intercept may vary between the subjects (random intercept). So different subjects can have a different base level of $FA$ but the effects that switching from right to left side of the brain must be the same for each subject. A simple example: Assume I want to estimate the average weight gain per year for members of my family, so I measure my son and my cat at birth (don't think to much about this, it's just an example for illustration), year 1, etc. So if I then use weight ~ age + (1|subject), this model captures the fact that my cat was much smaller than my son at the time of birth. However, this model does not capture the fact that my son also grew much faster, so it will try to fit the same growth rate for both my cat and my son.
So you should also include random intercepts. Random intercepts can only be used on within-subject variables (think about this for a while, if it isn't obvious why this is the case than re-read the previous until you understand what random-intercept mean).
So you could use a model like FA ~ Age + Wholebrain_FA + Volume + Group + Tract + Hemisphere + (Tract + Hemisphere|Subjects). In this model, the effects of Tract and Hemisphere are also allowed to vary between subjects. If I understand your data correctly, the tracts have been measured for both hemispheres. So you could even nest factors: FA ~ Age + Wholebrain_FA + Volume + Group + Tract + Hemisphere + (Tract|Hemisphere:Subjects) + (Hemisphere | Subjects), which indicates that the effect of tract could vary for each hemisphere within each subject, and the effect of hemisphere can vary between subjects. 
Warning: If you don't have enough data, a complicated model often cannot be computed. If you include all random intercepts and random slopes, the model may not converge. Nevertheless, it is often a good idea to start with the most complete model (all random intercepts and slopes) and then work from there to reduce parts of the model until you get convergence. Have a look at the between-subject variances and correlations of between-subject differences to see what should be dropped (often it is advisable to drop the smallest between-subject effects first) until you reach convergence. You should also describe this process in your publication: "We started with a model that included all random intercepts and slopes. Because this model did not converge, we dropped the random slope on Factor XYZ (SD=...). After dropping this random slope, the final model converged."
I found this page to be a very good primer on multilevel models in R.
