Quantify strength of association of two continuous variables while controlling for random effects I have a data set from a repeated measures experimental design with different sets of stimuli. I want to know how strong the association between the continuous dependent variable and the continuous predictor is while accounting for the interindividual and interstimulus variation.
My lmer model description in R looks like this
dv ~ pred + (1 | subject) + (1 | stimulus) 

Question 1: I understand that it is non-trivial to calculate R squared for random intercept-slope models. Is the same true for random intercept models? Is there an R-implementation of any of the available methods?
Question 2: If I z-transform my dependent variable and predictor, will the parameter estimate for the fixed effect reflect the strength of the association such as it would in an ordinary regression model? Update: I think I phrased this questions too vaguely. I was wondering if scaling variables would yield standardized regression estimates. I found this question has been answered before in another question.
Question 3: Is there an entirely different/more appropriate way to quantify the association strength while controlling for the interindividual and interstimulus variation?
 A: Q1: Yes, it is not trivial.  Please check the following link http://glmm.wikidot.com/faq under the section How do I compute a coefficient of determination (R2), or an analogue, for (G)LMMs?  as provides a quite good treatment of the issue; small 3-4 liners of code are presented that calculate an $R^2$ like measure. 
A rather helpful reference of the subject seems to be Xu's: Measuring explained variation in linear mixed effects models (2003) where the author computes the residual variance of the full model against the residual variance of a (fixed) intercept-only null model (in R syntax:  1-var(residuals(m))/(var(model.response(model.frame(m)))).
Q2: No, it would not reflect that. A crossed random effects structure in your data will still apply to the error structure you would assume in an OLS; you will just change the variances in your sample. 
Just to put things into perspectice: your OLS assumes that  $y \sim N(X\beta,\sigma^2 I)$, but you actually have that error diagonal structure only if you account for the "grouping in your samples" through the indicator variables $\gamma$, ie. $y|\gamma \sim N(X\beta+Z\gamma,\sigma^2I)$ (or unconditionally   $y\sim N(X\beta,\sigma^2I+ZDZ^T)$), whitening your $X$ and $y$ does change the amplitude of $D$ and $\sigma$ ($D$ is the diagonal matrix holding the std. deviation of your random effects) but other than that your structure remains unchanged.
Unless you have problems with numerical stability I would not recommend z-transforming / whitening your sample variables. 
Q3: Some people might suggest Generalized Estimating Equations (GEEs) but I think your approach is rather sensible and straightforward. 
