# wolfram-mathematica – DNA序列的混沌游戏

http://facstaff.unca.edu/mcmcclur/blog/GeneCGR.html

```genome = Import["c:\data\sequence.fasta", "Sequence"];
genome = StringReplace[ToString[genome], {"{" -> "", "}" -> ""}];
chars = StringCases[genome, "G" | "C" | "T" | "A"];
f[x_, "A"] := x/2;
f[x_, "T"] := x/2 + {1/2, 0};
f[x_, "G"] := x/2 + {1/2, 1/2};
f[x_, "C"] := x/2 + {0, 1/2};
pts = FoldList[f, {0.5, 0.5}, chars];
Graphics[{PointSize[Tiny], Point[pts]}]```

```shifts = chars /. {"A" -> {0., 0.}, "T" -> {.5, 0.}, "G" -> {.5, .5}, "C" -> {0, .5}};
fun1d = Compile[{{a, _Real, 1}}, FoldList[#/2 + #2 &, .5, a], CompilationTarget -> "C"]
pts = Transpose[fun1d /@ Transpose[shifts]];```

```threshold = 1;
With[{size = 300},
Image[1 - UnitStep[BinCounts[pts, 1/size, 1/size] - threshold]]
]```

```chars = RandomChoice[{"A", "T", "C", "G"}, 800000];

f[x_, "A"] := x/2;
f[x_, "T"] := x/2 + {1/2, 0};
f[x_, "G"] := x/2 + {1/2, 1/2};
f[x_, "C"] := x/2 + {0, 1/2};
Timing[pts = FoldList[f, {0.5, 0.5}, chars];]
Graphics[{PointSize[Tiny], Point[pts]}]```

```With[{resolution = 0.01},
ArrayPlot@BinCounts[pts, resolution, resolution]
]```

`arr = chars /. {"A" -> {0., 0.}, "T" -> {.5, 0.}, "G" -> {.5, .5}, "C" -> {0, .5}};`

```fun = Compile[{{a, _Real, 2}}, FoldList[#/2 + #2 &, {.5, .5}, a],
CompilationTarget -> "C"]```

`fun[arr]; // Timing`

```fun1d = Compile[{{a, _Real, 1}}, FoldList[#/2 + #2 &, .5, a],
CompilationTarget -> "C"]

arrt = Transpose[arr];
Timing[result = fun1d /@ arrt;]
pts = Transpose[result];```