Authors: J. Coelho, S. Ammer

library(ggplot2)
library(Momocs) # package for outline and curve analysis
This is Momocs 1.1.0

Attaching package: ‘Momocs’

The following object is masked from ‘package:stats’:

    filter

The following object is masked from ‘package:base’:

    table
setwd("/Users/del/ACADEMY/WIP/HumerusProject/Photos_final/") # location of my data
rawLPF <- import_tps(tps.path = "LPF/LPF.TPS", curves = T) # read data Left, Posterior, Females
rawLPM <- import_tps(tps.path = "LPM/LPM.TPS", curves = T) # read data Left, Posterior, Males

Rebuild dataset

For some reason, I can’t use Momocs::Out() directly on my imported TPS. Tried many ways, always obtaining different errors. I thought this might be because the p of coordinates is not consistent. So I resample them. In my understanding coo_sample() works well for open curves, while I preferred coo_interpolate() for my closed outlines, since they had very inconsistent p.

# Rebuild datasets
# For females
for (i in seq_along(rawLPF$cur)){
    # open (curves)
    rawLPF$cur[[i]][[1]] <- coo_sample(rawLPF$cur[[i]][[1]], n = 16)
    # closed (outlines)
    rawLPF$cur[[i]][[2]] <- coo_interpolate(rawLPF$cur[[i]][[2]], n = 24)
    # keep the IDs
    names(rawLPF$cur)[[i]] <- names(rawLPF$cur)[[i]]
}
## For males
for (i in seq_along(rawLPM$cur)){
    # open (curves)
    rawLPM$cur[[i]][[1]] <- coo_sample(rawLPM$cur[[i]][[1]], n = 16)
    # closed (outlines)
    rawLPM$cur[[i]][[2]] <- coo_interpolate(rawLPM$cur[[i]][[2]], n = 24)
    # keep the IDs
    names(rawLPM$cur)[[i]] <- names(rawLPM$cur)[[i]]
}

Spliting data

I tried to work with my curves together but for some reason, I also got errors trying to Out() these. So first I will split my open curve (TE: Throclear Extension) from my closed outlines (OF: Olecranon Fossa), but combine my males and females datasets to see the differences among these groups.

TE.f <- list()
for (i in seq_along(rawLPF$cur)){
    TE.f[[i]] <- rawLPF$cur[[i]][[1]]
    names(TE.f)[[i]] <- names(rawLPF$cur)[[i]]
}
TE.m <- list()
for (i in seq_along(rawLPM$cur)){
    TE.m[[i]] <- rawLPM$cur[[i]][[1]]
    names(TE.m)[[i]] <- names(rawLPM$cur)[[i]]
}
Sex <- rep("Female", length(TE.f))
TE.f <- Out(TE.f, fac = as.data.frame(Sex))
Sex <- rep("Male", length(TE.m))
TE.m <- Out(TE.m, fac = as.data.frame(Sex))
TE <- combine(TE.f, TE.m)
###
OF.f <- list()
for (i in seq_along(rawLPF$cur)){
    OF.f[[i]] <- rawLPF$cur[[i]][[2]]
    names(OF.f)[[i]] <- names(rawLPF$cur)[[i]]
}
OF.m <- list()
for (i in seq_along(rawLPM$cur)){
    OF.m[[i]] <- rawLPM$cur[[i]][[2]]
    names(OF.m)[[i]] <- names(rawLPM$cur)[[i]]
}
Sex <- rep("Female", length(OF.f))
OF.f <- Out(OF.f, fac = as.data.frame(Sex))
Sex <- rep("Male", length(OF.m))
OF.m <- Out(OF.m, fac = as.data.frame(Sex))
OF <- combine(OF.f, OF.m)

Procrustes and Rotation

My pics were (anatomically) upside-down. So I use the very useful coo_rotate() by Pi radians (180º). I have done some landmarks-based analysis before were Procrustes was a mandatory step. I am not sure if it is with open curves and closed outlines, however I am also doing a superimposition step here.

TE.aligned <- TE %>% coo_rotate(.,theta = pi) %>% fgProcrustes() %>% Opn()
OF.aligned <- OF %>% coo_rotate(.,theta = pi) %>% fgProcrustes() %>% coo_close()

Throclear Extension

MANOVA_PW(TE.pca, 1)
PC axes 1 to 2 were retained
$stars.tab
       Female Male
Female        .   

$summary (see also $manovas)
              Df  Pillai approx F num Df den Df Pr(>F)
Female - Male  1 0.03713    2.854      2    148 0.0608

Olecranon Fossa

OF.lda <- LDA(OF.pca, 1)
0.99 total variance
7 PC retained
LS0tCnRpdGxlOiAiRGlzdGFsIEh1bWVyaSBTaGFwZSBBbmFseXNpcyBOb3RlYm9vayIKb3V0cHV0OgogIHdvcmRfZG9jdW1lbnQ6IGRlZmF1bHQKICBodG1sX25vdGVib29rOiBkZWZhdWx0CiAgaHRtbF9kb2N1bWVudDogZGVmYXVsdAotLS0KCipBdXRob3JzOiogSi4gQ29lbGhvLCBTLiBBbW1lcgoKYGBge3IgaW1wb3J0aW5nIGRhdGEsIGNhY2hlPVRSVUV9CmxpYnJhcnkoZ2dwbG90MikKbGlicmFyeShNb21vY3MpICMgcGFja2FnZSBmb3Igb3V0bGluZSBhbmQgY3VydmUgYW5hbHlzaXMKc2V0d2QoIi9Vc2Vycy9kZWwvQUNBREVNWS9XSVAvSHVtZXJ1c1Byb2plY3QvUGhvdG9zX2ZpbmFsLyIpICMgbG9jYXRpb24gb2YgbXkgZGF0YQoKcmF3TFBGIDwtIGltcG9ydF90cHModHBzLnBhdGggPSAiTFBGL0xQRi5UUFMiLCBjdXJ2ZXMgPSBUKSAjIHJlYWQgZGF0YSBMZWZ0LCBQb3N0ZXJpb3IsIEZlbWFsZXMKcmF3TFBNIDwtIGltcG9ydF90cHModHBzLnBhdGggPSAiTFBNL0xQTS5UUFMiLCBjdXJ2ZXMgPSBUKSAjIHJlYWQgZGF0YSBMZWZ0LCBQb3N0ZXJpb3IsIE1hbGVzCgpgYGAKCiMjIFJlYnVpbGQgZGF0YXNldAoKRm9yIHNvbWUgcmVhc29uLCBJIGNhbid0IHVzZSBgTW9tb2NzOjpPdXQoKWAgZGlyZWN0bHkgb24gbXkgaW1wb3J0ZWQgVFBTLiBUcmllZCBtYW55IHdheXMsIGFsd2F5cyBvYnRhaW5pbmcgZGlmZmVyZW50IGVycm9ycy4gSSB0aG91Z2h0IHRoaXMgbWlnaHQgYmUgYmVjYXVzZSB0aGUgKipwKiogb2YgY29vcmRpbmF0ZXMgaXMgbm90IGNvbnNpc3RlbnQuIFNvIEkgcmVzYW1wbGUgdGhlbS4gSW4gbXkgdW5kZXJzdGFuZGluZyBgY29vX3NhbXBsZSgpYCB3b3JrcyB3ZWxsIGZvciBvcGVuIGN1cnZlcywgd2hpbGUgSSBwcmVmZXJyZWQgYGNvb19pbnRlcnBvbGF0ZSgpYCBmb3IgbXkgY2xvc2VkIG91dGxpbmVzLCBzaW5jZSB0aGV5IGhhZCB2ZXJ5IGluY29uc2lzdGVudCAqKnAqKi4KCmBgYHtyIHJlYnVpbGQsIGNhY2hlPVRSVUV9CiMgUmVidWlsZCBkYXRhc2V0cwojIEZvciBmZW1hbGVzCgpmb3IJKGkgaW4gc2VxX2Fsb25nKHJhd0xQRiRjdXIpKXsKCSMgb3BlbiAoY3VydmVzKQoJcmF3TFBGJGN1cltbaV1dW1sxXV0gPC0gY29vX3NhbXBsZShyYXdMUEYkY3VyW1tpXV1bWzFdXSwgbiA9IDE2KQoJIyBjbG9zZWQgKG91dGxpbmVzKQoJcmF3TFBGJGN1cltbaV1dW1syXV0gPC0gY29vX2ludGVycG9sYXRlKHJhd0xQRiRjdXJbW2ldXVtbMl1dLCBuID0gMjQpCgkjIGtlZXAgdGhlIElEcwoJbmFtZXMocmF3TFBGJGN1cilbW2ldXSA8LSBuYW1lcyhyYXdMUEYkY3VyKVtbaV1dCn0KCiMjIEZvciBtYWxlcwoKZm9yCShpIGluIHNlcV9hbG9uZyhyYXdMUE0kY3VyKSl7CgkjIG9wZW4gKGN1cnZlcykKCXJhd0xQTSRjdXJbW2ldXVtbMV1dIDwtIGNvb19zYW1wbGUocmF3TFBNJGN1cltbaV1dW1sxXV0sIG4gPSAxNikKCSMgY2xvc2VkIChvdXRsaW5lcykKCXJhd0xQTSRjdXJbW2ldXVtbMl1dIDwtIGNvb19pbnRlcnBvbGF0ZShyYXdMUE0kY3VyW1tpXV1bWzJdXSwgbiA9IDI0KQoJIyBrZWVwIHRoZSBJRHMKCW5hbWVzKHJhd0xQTSRjdXIpW1tpXV0gPC0gbmFtZXMocmF3TFBNJGN1cilbW2ldXQp9CgpgYGAKCiMjIFNwbGl0aW5nIGRhdGEKCkkgdHJpZWQgdG8gd29yayB3aXRoIG15IGN1cnZlcyB0b2dldGhlciBidXQgZm9yIHNvbWUgcmVhc29uLCBJIGFsc28gZ290IGVycm9ycyB0cnlpbmcgdG8gYE91dCgpYCB0aGVzZS4gU28gZmlyc3QgSSB3aWxsIHNwbGl0IG15IG9wZW4gY3VydmUgKFRFOiBUaHJvY2xlYXIgRXh0ZW5zaW9uKSBmcm9tIG15IGNsb3NlZCBvdXRsaW5lcyAoT0Y6IE9sZWNyYW5vbiBGb3NzYSksIGJ1dCBjb21iaW5lIG15IG1hbGVzIGFuZCBmZW1hbGVzIGRhdGFzZXRzIHRvIHNlZSB0aGUgZGlmZmVyZW5jZXMgYW1vbmcgdGhlc2UgZ3JvdXBzLgoKYGBge3Igc3BsaXQsIGNhY2hlPVRSVUV9CgpURS5mIDwtIGxpc3QoKQpmb3IJKGkgaW4gc2VxX2Fsb25nKHJhd0xQRiRjdXIpKXsKCVRFLmZbW2ldXSA8LSByYXdMUEYkY3VyW1tpXV1bWzFdXQoJbmFtZXMoVEUuZilbW2ldXSA8LSBuYW1lcyhyYXdMUEYkY3VyKVtbaV1dCn0KClRFLm0gPC0gbGlzdCgpCmZvcgkoaSBpbiBzZXFfYWxvbmcocmF3TFBNJGN1cikpewoJVEUubVtbaV1dIDwtIHJhd0xQTSRjdXJbW2ldXVtbMV1dCgluYW1lcyhURS5tKVtbaV1dIDwtIG5hbWVzKHJhd0xQTSRjdXIpW1tpXV0KfQoKU2V4IDwtIHJlcCgiRmVtYWxlIiwgbGVuZ3RoKFRFLmYpKQpURS5mIDwtIE91dChURS5mLCBmYWMgPSBhcy5kYXRhLmZyYW1lKFNleCkpClNleCA8LSByZXAoIk1hbGUiLCBsZW5ndGgoVEUubSkpClRFLm0gPC0gT3V0KFRFLm0sIGZhYyA9IGFzLmRhdGEuZnJhbWUoU2V4KSkKVEUgPC0gY29tYmluZShURS5mLCBURS5tKQoKIyMjCgpPRi5mIDwtIGxpc3QoKQpmb3IJKGkgaW4gc2VxX2Fsb25nKHJhd0xQRiRjdXIpKXsKCU9GLmZbW2ldXSA8LSByYXdMUEYkY3VyW1tpXV1bWzJdXQoJbmFtZXMoT0YuZilbW2ldXSA8LSBuYW1lcyhyYXdMUEYkY3VyKVtbaV1dCn0KCk9GLm0gPC0gbGlzdCgpCmZvcgkoaSBpbiBzZXFfYWxvbmcocmF3TFBNJGN1cikpewoJT0YubVtbaV1dIDwtIHJhd0xQTSRjdXJbW2ldXVtbMl1dCgluYW1lcyhPRi5tKVtbaV1dIDwtIG5hbWVzKHJhd0xQTSRjdXIpW1tpXV0KfQoKU2V4IDwtIHJlcCgiRmVtYWxlIiwgbGVuZ3RoKE9GLmYpKQpPRi5mIDwtIE91dChPRi5mLCBmYWMgPSBhcy5kYXRhLmZyYW1lKFNleCkpClNleCA8LSByZXAoIk1hbGUiLCBsZW5ndGgoT0YubSkpCk9GLm0gPC0gT3V0KE9GLm0sIGZhYyA9IGFzLmRhdGEuZnJhbWUoU2V4KSkKT0YgPC0gY29tYmluZShPRi5mLCBPRi5tKQoKCmBgYAoKIyMgUHJvY3J1c3RlcyBhbmQgUm90YXRpb24KCk15IHBpY3Mgd2VyZSAoYW5hdG9taWNhbGx5KSB1cHNpZGUtZG93bi4gU28gSSB1c2UgdGhlIHZlcnkgdXNlZnVsIGBjb29fcm90YXRlKClgIGJ5IFBpIHJhZGlhbnMgKDE4MMK6KS4gSSBoYXZlIGRvbmUgc29tZSBsYW5kbWFya3MtYmFzZWQgYW5hbHlzaXMgYmVmb3JlIHdlcmUgUHJvY3J1c3RlcyB3YXMgYSBtYW5kYXRvcnkgc3RlcC4gSSBhbSBub3Qgc3VyZSBpZiBpdCBpcyB3aXRoIG9wZW4gY3VydmVzIGFuZCBjbG9zZWQgb3V0bGluZXMsIGhvd2V2ZXIgSSBhbSBhbHNvIGRvaW5nIGEgc3VwZXJpbXBvc2l0aW9uIHN0ZXAgaGVyZS4KCmBgYHtyIFByb2NydXN0ZXMsIGNhY2hlPVRSVUV9ClRFLmFsaWduZWQgPC0gVEUgJT4lIGNvb19yb3RhdGUoLix0aGV0YSA9IHBpKSAlPiUgZmdQcm9jcnVzdGVzKCkgJT4lIE9wbigpCk9GLmFsaWduZWQgPC0gT0YgJT4lIGNvb19yb3RhdGUoLix0aGV0YSA9IHBpKSAlPiUgZmdQcm9jcnVzdGVzKCkgJT4lIGNvb19jbG9zZSgpCmBgYAoKIyBUaHJvY2xlYXIgRXh0ZW5zaW9uCgpgYGB7ciBURSBwbG90dGluZywgY2FjaGU9VFJVRSwgZmlnLndpZHRoPTEyfQoKVEUuYWxpZ25lZCAlPiUgc3RhY2soLiwgdGl0bGUgPSAiVHJvY2hsZWFyIEV4dGVuc2lvbiAoSC4gc2FwaWVucykiLCBmYWMgPSAiU2V4IikKZnRpbWVzIDwtIHJlcCgncmVkJywgbGVuZ3RoKGdyZXAoJ0ZlbWFsZScsIFRFLmFsaWduZWQkZmFjJFNleCkpKQptdGltZXMgPC0gcmVwKCdibHVlJywgbGVuZ3RoKGdyZXAoJ01hbGUnLCBURS5hbGlnbmVkJGZhYyRTZXgpKSkKcGFuZWwoVEUuYWxpZ25lZCwgYm9yZGVycyA9IGMoZnRpbWVzLCBtdGltZXMpLCBuYW1lcyA9IHN1YnN0cihuYW1lcyhURS5hbGlnbmVkKSwgNywgOSkpCgpURS5hbGlnbmVkICU+JSBucG9seSgpICU+JSBQQ0EoKSAlPiUgcGxvdCguLCJTZXgiKSAjIEkgbmVlZCB0byByZWFkIGFib3V0IHRoaXMKVEUuYWxpZ25lZCAlPiUgb3BvbHkoKSAlPiUgUENBKCkgJT4lIHBsb3QoLiwiU2V4IikgIyByZWFkIGFib3V0IHRoaXMKVEUuYWxpZ25lZCAlPiUgZGZvdXJpZXIoKSAlPiUgUENBKCkgJT4lIHBsb3QoLiwiU2V4IikgIyBzdHJhbmdlIHJlc3VsdHM7IG1pZ2h0IG5vdCBiZSBwcm9wZXIgbWV0aG9kIGZvciB0aGlzIGRhdGFzZXQuCgpURS5uIDwtIG5wb2x5KFRFLmFsaWduZWQsIG5iLnB0cyA9IDE2LCBkZWdyZWUgPSA1KQpURS5wY2EgPC0gUENBKFRFLm4pClRFLnBjYSAlPiUgYXNfZGYoKSAlPiUgZ2dwbG90KCkgKwoJYWVzKHggPSBQQzEsIHkgPSBQQzIsIGNvbCA9IFNleCkgKyBjb29yZF9lcXVhbCgpICsgCglnZW9tX3BvaW50KCkgKyBnZW9tX2RlbnNpdHkyZCgpICsgdGhlbWVfbGlnaHQoKQoKIyBtZWFuIHNoYXBlOgpURS5uICU+JSBtc2hhcGVzKCkgJT4lIGNvb19wbG90KCkKIyBtZWFuIHNoYXBlIHBlciBncm91cDoKVEUubXMgPC0gbXNoYXBlcyhURS5uLCAxKQpURV9mZW1hbGUgPC0gVEUubXMkc2hwJEZlbWFsZSAlVD4lIGNvb19wbG90KGJvcmRlciA9ICJyZWQiKQpURV9tYWxlIDwtIFRFLm1zJHNocCRNYWxlICVUPiUgY29vX2RyYXcoYm9yZGVyID0gImJsdWUiKQpsZWdlbmQoInRvcHJpZ2h0IiwgbHdkID0gMSwgY29sID0gYygicmVkIiwgImJsdWUiKSwgbGVnZW5kID0gYygiRmVtYWxlIiwgIk1hbGUiKSwgY2V4ID0gMC44KQp0aXRsZShtYWluID0gIlRocm9jbGVhciBFeHRlbnNpb24gLSBNZWFuIHNoYXBlcyIpCgojIFRQUyBwbG90cwpjb29fYXJyb3dzKFRFX2ZlbWFsZSwgVEVfbWFsZSk7IHRpdGxlKCJEZWZvcm1hdGlvbnMiKQp0cHNfZ3JpZChURV9mZW1hbGUsIFRFX21hbGUpCnRwc19hcnIoVEVfZmVtYWxlLCBURV9tYWxlKQp0cHNfaXNvKFRFX2ZlbWFsZSwgVEVfbWFsZSkKCiMgTWFub3ZhIGFuYWx5c2lzCgpNQU5PVkEoVEUucGNhLCAxKQpNQU5PVkFfUFcoVEUucGNhLCAxKQoKIyBMREEgLyBDVkEgYW5hbHlzaXMKClRFLmxkYSA8LSBMREEoVEUucGNhLCAxKQpURS5sZGEKcGxvdChURS5sZGEpCnBsb3RfQ1YoVEUubGRhKQpwbG90X0NWMihURS5sZGEpCgpgYGAKCiMgT2xlY3Jhbm9uIEZvc3NhCgpgYGB7ciBPRiBwbG90dGluZywgY2FjaGU9VFJVRSwgZmlnLndpZHRoPTEyfQoKT0YuYWxpZ25lZCAlPiUgc3RhY2soLix0aXRsZSA9ICJPbGVjcmFub24gRm9zc2EgKEguIHNhcGllbnMpIiwgZmFjID0gIlNleCIpCk9GLmFsaWduZWQgJT4lIHBhbmVsKC4sIGZhYyA9ICJTZXgiLCBuYW1lcyA9IHN1YnN0cihuYW1lcyhPRi5hbGlnbmVkKSwgNywgOSkpCgpPRi5hbGlnbmVkICU+JSBlZm91cmllcihub3JtPVRSVUUpICU+JSBQQ0EoKSAlPiUgcGxvdCguLCJTZXgiKQoKT0YuZWYgPC0gZWZvdXJpZXIoT0YuYWxpZ25lZCwgbmIuaCA9IDYpCk9GLnBjYSA8LSBQQ0EoT0YuZWYpCk9GLnBjYSAlPiUgYXNfZGYoKSAlPiUgZ2dwbG90KCkgKwoJYWVzKHggPSBQQzEsIHkgPSBQQzIsIGNvbCA9IFNleCkgKyBjb29yZF9lcXVhbCgpICsgCglnZW9tX3BvaW50KCkgKyBnZW9tX2RlbnNpdHkyZCgpICsgdGhlbWVfbGlnaHQoKQoKIyBtZWFuIHNoYXBlOgpPRi5lZiAlPiUgbXNoYXBlcygpICU+JSBjb29fcGxvdCgpCiMgbWVhbiBzaGFwZSBwZXIgZ3JvdXA6Ck9GLm1zIDwtIG1zaGFwZXMoT0YuZWYsIDEpCk9GX2ZlbWFsZSA8LSBPRi5tcyRzaHAkRmVtYWxlICVUPiUgY29vX3Bsb3QoYm9yZGVyID0gInJlZCIpCk9GX21hbGUgPC0gT0YubXMkc2hwJE1hbGUgJVQ+JSBjb29fZHJhdyhib3JkZXIgPSAiYmx1ZSIpCmxlZ2VuZCgidG9wcmlnaHQiLCBsd2QgPSAxLCBjb2wgPSBjKCJyZWQiLCAiYmx1ZSIpLCBsZWdlbmQgPSBjKCJGZW1hbGUiLCAiTWFsZSIpLCBjZXggPSAwLjgpCnRpdGxlKG1haW4gPSAiT2xlY3Jhbm9uIEZvc3NhIC0gTWVhbiBzaGFwZXMiKQoKIyBUUFMgcGxvdHMKY29vX2Fycm93cyhPRl9mZW1hbGUsIE9GX21hbGUpOyB0aXRsZSgiRGVmb3JtYXRpb25zIikKdHBzX2dyaWQoT0ZfZmVtYWxlLCBPRl9tYWxlKQp0cHNfYXJyKE9GX2ZlbWFsZSwgT0ZfbWFsZSkKdHBzX2lzbyhPRl9mZW1hbGUsIE9GX21hbGUsIGdyaWQgPSBUUlVFKQoKIyBNQU5PVkEgYW5hbHlzaXMKCk1BTk9WQShPRi5wY2EsIDEpCk1BTk9WQV9QVyhPRi5wY2EsIDEpCgojIExEQSAvIENWQSBhbmFseXNpcwoKT0YubGRhIDwtIExEQShPRi5wY2EsIDEpCk9GLmxkYQpwbG90KE9GLmxkYSkKcGxvdF9DVihPRi5sZGEpCnBsb3RfQ1YyKE9GLmxkYSkKCgoKYGBgCgo=