]> git.vanrenterghem.biz Git - www2.vanrenterghem.biz.git/blob - source/posts/Fibonacci_golden_spiral.org
90abf9a3f439344d5e9f7b94d1c46db6d3e4e4d1
[www2.vanrenterghem.biz.git] / source / posts / Fibonacci_golden_spiral.org
1 #+TITLE: Creating a golden spiral in R\r
2 #+AUTHOR: Frederik Vanrenterghem\r
3 #+filetags: :R:code:\r
4 #+LANGUAGE: en\r
5 #+PROPERTY: session *R* \r
6 #+PROPERTY: cache yes \r
7 #+PROPERTY: results file graphics \r
8 #+PROPERTY: exports both \r
9 #+PROPERTY: tangle yes \r
10 #+OPTIONS: toc:nil  \r
11 #+date: <2019-09-16 22:03:03>\r
12 \r
13 * What \r
14 After having read the first part of a Rcpp tutorial which compared native R vs C++ implementations of a Fibonacci sequence generator, I resorted to drawing the so-called Golden Spiral using R.\r
15 \r
16 * Details\r
17 Libraries used in this example are the following\r
18 \r
19 #+BEGIN_SRC R :session R\r
20 library(ggplot2)\r
21 library(plotrix)\r
22 #+END_SRC\r
23 \r
24 #+RESULTS:\r
25 | plotrix   |\r
26 | ggplot2   |\r
27 | stats     |\r
28 | graphics  |\r
29 | grDevices |\r
30 | utils     |\r
31 | datasets  |\r
32 | methods   |\r
33 | base      |\r
34 \r
35 In polar coordinates, this special instance of a logarithmic spiral's functional representation can be simplified to r(t) = e^(0.0635*t)\r
36 For every quarter turn, the corresponding point on the spiral is a factor of phi further from the origin (r is this distance), with phi the golden ratio - the same one obtained from dividing any 2 sufficiently big successive numbers on a Fibonacci sequence, which is how the golden ratio, the golden spiral and Fibonacci sequences are linked concepts!\r
37 #+BEGIN_SRC R :session R\r
38 polar_golden_spiral <- function(theta) exp(0.30635*theta)\r
39 #+END_SRC\r
40 \r
41 #+RESULTS:\r
42 \r
43 Let's do 2 full circles. First, I create a sequence of angle values theta. Since 2 * PI is the equivalent of a circle in polar coordinates, we need to have distances from origin for values between 0 and 4 * PI.\r
44 #+BEGIN_SRC R :session R\r
45 seq_theta <- seq(0,4*pi,by=0.05)\r
46 \r
47 dist_from_origin <- sapply(seq_theta,polar_golden_spiral)\r
48 #+END_SRC\r
49 \r
50 #+RESULTS:\r
51 |                1 |\r
52 | 1.01543541418402 |\r
53 | 1.03110908037907 |\r
54 | 1.04702467610363 |\r
55 | 1.06318593564018 |\r
56 | 1.07959665091141 |\r
57 | 1.09626067236991 |\r
58 | 1.11318190990159 |\r
59 | 1.13036433374308 |\r
60 | 1.14781197541325 |\r
61 | 1.16552892865913 |\r
62 | 1.18351935041644 |\r
63 | 1.20178746178492 |\r
64 | 1.22033754901874 |\r
65 | 1.23917396453215 |\r
66 | 1.25830112792076 |\r
67 | 1.27772352699844 |\r
68 | 1.29744571885033 |\r
69 | 1.31747233090207 |\r
70 | 1.33780806200553 |\r
71 | 1.35845768354131 |\r
72 | 1.37942604053823 |\r
73 | 1.40071805281016 |\r
74 | 1.42233871611032 |\r
75 | 1.44429310330345 |\r
76 | 1.46658636555607 |\r
77 | 1.48922373354506 |\r
78 |   1.512210518685 |\r
79 | 1.53555211437434 |\r
80 | 1.55925399726085 |\r
81 | 1.58332172852666 |\r
82 | 1.60776095519303 |\r
83 | 1.63257741144533 |\r
84 | 1.65777691997847 |\r
85 | 1.68336539336305 |\r
86 | 1.70934883543265 |\r
87 | 1.73573334269253 |\r
88 |    1.76252510575 |\r
89 | 1.78973041076699 |\r
90 | 1.81735564093491 |\r
91 | 1.84540727797241 |\r
92 | 1.87389190364612 |\r
93 | 1.90281620131498 |\r
94 | 1.93218695749834 |\r
95 | 1.96201106346829 |\r
96 | 1.99229551686655 |\r
97 | 2.02304742334636 |\r
98 | 2.05427399823962 |\r
99 | 2.08598256824992 |\r
100 |  2.1181805731715 |\r
101 | 2.15087556763495 |\r
102 | 2.18407522287968 |\r
103 | 2.21778732855389 |\r
104 | 2.25201979454219 |\r
105 | 2.28678065282156 |\r
106 | 2.32207805934587 |\r
107 |  2.3579202959595 |\r
108 | 2.39431577234054 |\r
109 | 2.43127302797395 |\r
110 | 2.46880073415517 |\r
111 | 2.50690769602466 |\r
112 | 2.54560285463391 |\r
113 | 2.58489528904321 |\r
114 | 2.62479421845192 |\r
115 | 2.66530900436155 |\r
116 | 2.70644915277227 |\r
117 |  2.7482243164133 |\r
118 | 2.79064429700773 |\r
119 | 2.83371904757232 |\r
120 | 2.87745867475275 |\r
121 | 2.92187344119496 |\r
122 |  2.9669737679531 |\r
123 | 3.01277023693458 |\r
124 | 3.05927359338295 |\r
125 | 3.10649474839905 |\r
126 | 3.15444478150108 |\r
127 | 3.20313494322417 |\r
128 | 3.25257665776014 |\r
129 | 3.30278152563795 |\r
130 |  3.3537613264455 |\r
131 | 3.40552802159354 |\r
132 | 3.45809375712212 |\r
133 | 3.51147086655048 |\r
134 | 3.56567187377081 |\r
135 | 3.62070949598677 |\r
136 | 3.67659664669734 |\r
137 |  3.7333464387267 |\r
138 | 3.79097218730088 |\r
139 | 3.84948741317197 |\r
140 | 3.90890584579046 |\r
141 | 3.96924142652657 |\r
142 | 4.03050831194138 |\r
143 | 4.09272087710833 |\r
144 | 4.15589371898609 |\r
145 | 4.22004165984341 |\r
146 | 4.28517975073691 |\r
147 | 4.35132327504251 |\r
148 | 4.41848775204137 |\r
149 | 4.48668894056115 |\r
150 | 4.55594284267357 |\r
151 | 4.62626570744896 |\r
152 | 4.69767403476877 |\r
153 | 4.77018457919694 |\r
154 | 4.84381435391107 |\r
155 | 4.91858063469419 |\r
156 |  4.9945009639882 |\r
157 | 5.07159315500985 |\r
158 | 5.14987529593027 |\r
159 | 5.22936575411901 |\r
160 | 5.31008318045357 |\r
161 | 5.39204651369547 |\r
162 | 5.47527498493387 |\r
163 | 5.55978812209773 |\r
164 |  5.6456057545377 |\r
165 | 5.73274801767868 |\r
166 | 5.82123535774417 |\r
167 | 5.91108853655362 |\r
168 | 6.00232863639374 |\r
169 | 6.09497706496509 |\r
170 | 6.18905556040493 |\r
171 | 6.28458619638769 |\r
172 | 6.38159138730412 |\r
173 | 6.48009389352033 |\r
174 | 6.58011682671816 |\r
175 |  6.6816836553178 |\r
176 | 6.78481820998423 |\r
177 | 6.88954468921862 |\r
178 | 6.99588766503603 |\r
179 | 7.10387208873074 |\r
180 |  7.2135232967306 |\r
181 | 7.32486701654172 |\r
182 | 7.43792937278492 |\r
183 | 7.55273689332534 |\r
184 | 7.66931651549675 |\r
185 | 7.78769559242179 |\r
186 | 7.90790189942989 |\r
187 |  8.0299636405742 |\r
188 | 8.15390945524908 |\r
189 | 8.27976842490986 |\r
190 | 8.40757007989611 |\r
191 | 8.53734440636049 |\r
192 |  8.6691218533043 |\r
193 | 8.80293333972179 |\r
194 | 8.93881026185472 |\r
195 | 9.07678450055882 |\r
196 | 9.21688842878404 |\r
197 | 9.35915491917023 |\r
198 | 9.50361735176004 |\r
199 |  9.6503096218309 |\r
200 | 9.79926614784789 |\r
201 | 9.95052187953938 |\r
202 | 10.1041123060972 |\r
203 | 10.2600734645037 |\r
204 | 10.4184419479868 |\r
205 | 10.5792549146061 |\r
206 | 10.7425500959714 |\r
207 | 10.9083658060953 |\r
208 | 11.0767409503832 |\r
209 | 11.2477150347615 |\r
210 | 11.4213281749469 |\r
211 | 11.5976211058588 |\r
212 | 11.7766351911771 |\r
213 |  11.958412433047 |\r
214 | 12.1429954819344 |\r
215 | 12.3304276466328 |\r
216 | 12.5207529044246 |\r
217 | 12.7140159114002 |\r
218 | 12.9102620129349 |\r
219 | 13.1095372543288 |\r
220 | 13.3118883916102 |\r
221 | 13.5173629025061 |\r
222 |  13.726008997582 |\r
223 | 13.9378756315533 |\r
224 | 14.1530125147717 |\r
225 | 14.3714701248888 |\r
226 | 14.5932997186998 |\r
227 | 14.8185533441694 |\r
228 | 15.0472838526447 |\r
229 | 15.2795449112548 |\r
230 | 15.5153910155034 |\r
231 | 15.7548775020547 |\r
232 | 15.9980605617174 |\r
233 | 16.2449972526286 |\r
234 | 16.4957455136412 |\r
235 | 16.7503641779184 |\r
236 | 17.0089129867378 |\r
237 |  17.271452603508 |\r
238 | 17.5380446280028 |\r
239 | 17.8087516108139 |\r
240 | 18.0836370680272 |\r
241 | 18.3627654961257 |\r
242 | 18.6462023871224 |\r
243 | 18.9340142439267 |\r
244 | 19.2262685959479 |\r
245 | 19.5230340149396 |\r
246 | 19.8243801310889 |\r
247 | 20.1303776493537 |\r
248 | 20.4410983660522 |\r
249 | 20.7566151857085 |\r
250 | 21.0770021381583 |\r
251 | 21.4023343959182 |\r
252 | 21.7326882918241 |\r
253 | 22.0681413369407 |\r
254 | 22.4087722387478 |\r
255 | 22.7546609196083 |\r
256 | 23.1058885355194 |\r
257 |  23.462537495155 |\r
258 | 23.8246914792008 |\r
259 | 24.1924354599887 |\r
260 | 24.5658557214339 |\r
261 | 24.9450398792791 |\r
262 | 25.3300769016527 |\r
263 | 25.7210571299428 |\r
264 | 26.1180722999943 |\r
265 | 26.5212155636329 |\r
266 | 26.9305815105213 |\r
267 | 27.3462661903527 |\r
268 | 27.7683671353873 |\r
269 | 28.1969833833359 |\r
270 | 28.6322155005976 |\r
271 | 29.0741656058555 |\r
272 | 29.5229373940367 |\r
273 | 29.9786361606425 |\r
274 | 30.4413688264541 |\r
275 |  30.911243962619 |\r
276 | 31.3883718161253 |\r
277 | 31.8728643356692 |\r
278 | 32.3648351979214 |\r
279 | 32.8643998341988 |\r
280 |  33.371675457549 |\r
281 | 33.8867810902509 |\r
282 | 34.4098375917422 |\r
283 | 34.9409676869756 |\r
284 | 35.4802959952146 |\r
285 | 36.0279490592724 |\r
286 |  36.584055375203 |\r
287 | 37.1487454224504 |\r
288 | 37.7221516944627 |\r
289 | 38.3044087297792 |\r
290 | 38.8956531435973 |\r
291 | 39.4960236598267 |\r
292 |  40.105661143638 |\r
293 | 40.7247086345141 |\r
294 | 41.3533113798114 |\r
295 | 41.9916168688395 |\r
296 | 42.6397748674667 |\r
297 | 43.2979374532595 |\r
298 | 43.9662590511644 |\r
299 |  44.644896469741 |\r
300 | 45.3340089379542 |\r
301 | 46.0337581425336 |\r
302 | 46.7443082659106 |\r
303 \r
304 Plotting the function using coord_polar in ggplot2 does not work as intended. Unexpectedly, the x axis keeps extending instead of circling back once a full circle is reached. Turns out coord_polar might not really be intended to plot elements in polar vector format.\r
305 #+BEGIN_SRC R :session R :results output file graphics :file ../assets/golden_spiral-coord_polar-fail.png :exports both\r
306 ggplot(data.frame(x = seq_theta, y = dist_from_origin), aes(x,y)) +\r
307     geom_point() +\r
308     coord_polar(theta="x")\r
309 #+END_SRC\r
310 \r
311 #+RESULTS:\r
312 [[file:../assets/golden_spiral-coord_polar-fail.png]]\r
313 \r
314 To ensure what I was trying to do is possible, I employ a specialised plotfunction instead\r
315 #+BEGIN_SRC R :session R :results output file graphics :file ../assets/golden_spiral-plotrix.png :exports both\r
316 plotrix::radial.plot(dist_from_origin, seq_theta,rp.type="s", point.col = "blue")\r
317 #+END_SRC\r
318 \r
319 With that established and the original objective of the exercise achieved, it still would be nice to be able to accomplish this using ggplot2. To do so, the created sequence above needs to be converted to cartesian coordinates.\r
320 The rectangular function equivalent of the golden spiral function r(t) defined above is a(t) = (r(t) cos(t), r(t) sin(t))\r
321 It's not too hard to come up with a hack to convert one to the other.\r
322 #+BEGIN_SRC R :session R\r
323 cartesian_golden_spiral <- function(theta) {\r
324     a <- polar_golden_spiral(theta)*cos(theta)\r
325     b <- polar_golden_spiral(theta)*sin(theta)\r
326     c(a,b)\r
327 }\r
328 #+END_SRC\r
329 \r
330 #+RESULTS:\r
331 \r
332 Applying that function to the same series of angles from above and stitching the resulting coordinates in a data frame. Note I'm enclosing the first expression in brackets, which prints it immediately, which is useful when the script is run interactively.\r
333 #+BEGIN_SRC R :session R :exports code\r
334 (serie <- sapply(seq_theta,cartesian_golden_spiral))\r
335 df <- data.frame(t(serie))\r
336 #+END_SRC\r
337 \r
338 #+RESULTS:\r
339 : TRUE\r
340 \r
341 * Result\r
342 With everything now ready in the right coordinate system, it's now only a matter of setting some options to make the output look acceptable.\r
343 #+BEGIN_SRC R :session R :results output file graphics :file ../assets/golden_spiral-ggplot-coord-fixed.png :width 800 :height 800 :exports both\r
344 ggplot(df, aes(x=X1,y=X2)) +\r
345     geom_path(color="blue") +\r
346     theme(panel.grid.minor = element_blank(),\r
347           axis.text.x = element_blank(),\r
348           axis.text.y = element_blank()) +\r
349     scale_y_continuous(breaks = seq(-20,20,by=10)) +\r
350     scale_x_continuous(breaks = seq(-20,50,by=10)) +\r
351     coord_fixed() +\r
352     labs(title = "Golden spiral",\r
353          subtitle = "Another view on the Fibonacci sequence",\r
354          caption = "Maths from https://www.intmath.com/blog/mathematics/golden-spiral-6512\nCode errors mine.",\r
355          x = "",\r
356          y = "")\r
357 #+END_SRC\r
358 \r
359 * Note on how this post was written.\r
360 After a long hiatus, I set about using emacs, org-mode and ESS together to create this post. All code is part of an .org file, and gets exported to markdown using the orgmode conversion - C-c C-e m m.\r