Ali,
I think what you need is something like this (I modified an example from the documentation)
RSolve[{a[n + 1] - 2 a[n] + eps[n] == 1, a[0] == 1}, a[n], n]
Which returns a[n] as a function of n and eps[n]. Your eps[n] is a random sequence. To simulate it you would need to generate points for your eps using something like RandomVariate[].
I am not sure you will obtain solutions for complicated equations and you may just have to numerically generate those solutions using the recurrence relationship like this:
In[20]:= RecurrenceTable[{a[n + 1] == 3 a[n], a[1] == 7}, a, {n, 1,
10}]
Out[20]= {7, 21, 63, 189, 567, 1701, 5103, 15309, 45927, 137781}
In[22]:= RecurrenceTable[{a[n + 1] ==
3 a[n] + RandomVariate[NormalDistribution[1, 3]],
a[1] == 7}, a, {n, 1, 10}]
Out[22]= {7., 22.4549, 68.8196, 207.914, 625.196, 1877.04, 5632.58, \
16899.2, 50699., 152099.}