Unequal variance and crossed random effects in linear mixed effects model

Question

I'm analyzing an experiment that has 50 subjects and 50 items. Each item can occur in two possible conditions. Every experimental subject produces a (continuous) response to every item, but only once per item (individually, they only respond to one condition per item).

I hypothesize the following effects:

1. item condition should affect the mean response: higher response in one condition than the other
2. item condition should affect the variance of the response: higher variance in one condition than the other
3. there is idiosyncratic variation between subjects: some subjects produce a higher response than others, regardless of item or condition. In addition, some subjects are more sensitive to the item condition than others.
4. there is idiosyncratic variation between items: some items get higher responses than others, regardless of condition. In addition, some items are more sensitive to the item condition than others.

I want to model this with a linear mixed-effects model, and test the significance of (1) and (2). If I were using lme4, I would specify a model like this, with a fixed effect of condition, plus random intercepts and condition slopes for subjects and for items.

lmer(response ~ condition + (1 + condition | subject) + (1 + condition | item))

As far as I can tell, there's no way to test (2) with lme4. If someone can suggest a straightforward way to test (2) while also accounting for the other components, that would be very helpful!

This is what I've been able to do so far:

This question suggests how to use nlme to allow the variance to differ between fixed factor levels. I found via googling that I can specify crossed random effects with a hack using pdBlocked and pdIdent. My current model looks like this:

my_data$one = rep(1, nrow(my_data))

lme(response ~ condition, 
    random=list(one = pdBlocked(list(pdIdent(~ subject - 1),
                                     pdIdent(~ item - 1)))
                      ),
     weights=varIdent(form=~1|condition),
     data=my_data)

However, as far as I can tell from the model output (and comparison to lme4), this includes random intercepts but not slopes for condition, so requirements (3) and (4) of the model aren't fully satisfied. If you know how to specify condition slopes using the nlme syntax, that would also answer this question.