(* this function gives next position of the DNA after movement *)
NextPosition[{ROW_, COL_}, {x_Integer, y_Integer}] :=
  Module[{neighbors, length, distribution, weights, choice},
   (* get the neighbours, note that they are hard-
   coded in the order much that it matches the weights to be \
calculated later *)
   neighbors = {x, y} + # & /@ {{-1, 0}, {-1, 1}, {0, 1}, {1, 0}, {1, 
       1}};
   length = Length@neighbors;
   (* need to assign weight in such a way so that it favours the \
middle one with less emphasis on the extrem, 
   so binomial distribution with mean at center; note mean is n*p, 
   we know n *)
   distribution = 
    BinomialDistribution[length, Ceiling[length/2]/length];
   weights = (PDF[distribution, #] & /@ Range[length])*
     (* also need to check some of the neighbours are feasible *)
     (If[1 <= First@# <= ROW && 1 <= Last@# <= COL + 1, 1, 0] & /@ 
       neighbors);
   (* Print["weights\t"<>ToString[N@weights]]; *)
   choice = RandomChoice[weights -> neighbors];
   (* return Null to indicate left the region *)
   If[Last@choice < COL, choice, Null]];


(* this function is the one to be used in NestList with \
{positiveDNA_List, negativeDNA_List, tale_List}, also Sow successful \
bind *)
(* Transition[dna_List]:=Complement[NextPosition/@dna,tale]; *)
Transition[{ROW_, COL_}, {dnaPos_List, dnaNeg_List, tale_List, 
    taleBind_List}] := Module[{nextPos, nextNeg, binding},
   (* negative dna don't interact *)
   nextNeg = Select[NextPosition[{ROW, COL}, #] & /@ dnaNeg, ListQ];
   (* positive dna will always bind if in contact with tale *)
   nextPos = Select[NextPosition[{ROW, COL}, #] & /@ dnaPos, ListQ];
   binding = Intersection[nextPos, tale];
   (* Sow the count and return the result *)
   Sow[Length@binding];
   Return@{Complement[nextPos, binding], nextNeg, 
     Complement[tale, binding], Union[taleBind, binding]}
   ];

Simulation[{ROW_Integer, COL_Integer}, 
   time_Integer, {posP_, taleP_}] :=
  Module[{dnaPos, dnaNeg, tale, sim, result},
   (* initialize DNA position, 
   fill posP % of available position total *)
   dnaPos = 
    RandomSample[
     Table[{i, j}, {j, 1, Ceiling[0.45 COL]}, {i, 1, ROW}]~Flatten~1,
     Max[Floor[posP 0.45 0.5 ROW COL], 1]];
   dnaNeg = 
    RandomSample[
     Table[{i, j}, {j, 1, Ceiling[0.45 COL]}, {i, 1, ROW}]~Flatten~1,
     Max[Floor[(1 - posP) 0.45 0.5 ROW COL], 1]];
   (* initialize TALE position, fill 50% of available position *)
   tale = 
    RandomSample[
     Table[{i, j}, {i, 1, ROW}, {j, Floor[0.5 COL], COL}]~Flatten~1,
     Max[Floor[taleP 0.5 ROW COL], 1]];
   
   (* return the simulation result *)
   (* Return@Append[Transpose[{NestList[Transition,dnaPos,time],
   NestList[Transition,dnaNeg,time]}],tale]] *)
   result = 
    Last@Last@
      Reap[sim = 
        NestList[
         Transition[{ROW, COL}, #] &, {dnaPos, dnaNeg, tale, {}}, 
         time]];
   Return[{sim, result}]];

Manipulate[
 result = First@Simulation[{dim1, dim2}, time, {p1, p2}];
 Dynamic@MatrixPlot[
   SparseArray[
    Flatten@ (* negative dna (green), positive dna (red), 
     binding (blue), last is tale (black) *)
     {Thread[result[[n, 2]] -> 1],
      Thread[result[[n, 1]] -> 2],
      Thread[result[[n, 4]] -> 3],
      Thread[result[[n, 3]] -> 4]}, {dim1, dim2}],
   ColorRules -> {1 -> Green, 2 -> Red, 3 -> Blue, 4 -> Black}],
 "Global Settings",
 {{dim1, 20, "# of rows"}, 10, 40, 1},
 {{dim2, 20, "# of cols"}, 10, 40, 1},
 {{time, Max[dim1, dim2], "total # of steps"}, 10, 100, 1},
 {{p1, 0.5, "% filling of positive DNA"}, 0, 1},
 {{p2, 0.5, "% filling of TALEs"}, 0, 1},
 Delimiter,
 "Display settings",
 {{n, 1, "current step"}, 1, time, 1}, SaveDefinitions -> True, 
 ContinuousAction -> False, ControlPlacement -> Top]