function y=logistic(x,p,q)
% random logistic map. For this to work correctly, need 1>p>-=.5 and 0<q<1
y=[]; n=1:30000; x0=x;
for i=n
    x=x*(1-x*((1+q)*(rand<p)-q));
    %x=x*(1-x*((1+q)*p-q)); % deterministic analogue-don't use p=0.5
    y=[y x];
end
plot(log(n),log(y))
grid on
xlabel('ln(n)')
ylabel(['ln(x)    x_0=' num2str(x0)])
title(['Random logistic map X_{n+1}=X_n(1-X_n((1+q)I_n-q))   p= ' num2str(p) ' q=' num2str(q)])