jueves, 5 de mayo de 2011

Càlculo de distancias euclidianas a partir de landmarks en R.

NOTA: este post puede parecer muy largo, pero no lo es!!; lo que pasa es que se muestran los resultados para que sean comparados con los que obtendrás en casa al hacer la practica.

Muchas veces se necesitan realizar análisis de morfometrìa en los cuales aun se utilizan distancias euclidianas para llevar a cabo análisis de variaciones en tamaño como parte de estudios de tipo "tradicional".
Este tipo de análisis, a pesar de "boom" de el estudio de forma sin relación con el tamaño, esta aun vigente y es importante conocer aspectos básicos acerca del uso de técnicas de distancias "tradicionales".

Básicamente para el calculo de distancias euclidianas, se utiliza (valga la redundancia) la formula de distancias en el espacio euclidiano en un sistema de coordenadas .

Asi a partir de la toma de landmarks y los archivos que podamos tener de landmarks de varios individuos, podemos medir varias distancias en varios individuos con un solo Script, siguiendo una serie de comandos y ponerlos todos en un solo archivo. (los archivos todos deben haber sido tomados utilizando la misma escala (objeto de referencia de distancia conocida), puesto que estamos trabajando con tamaño)

De tal modo que empezare explicando el calculo de distancias en un solo individuo,(1) primero tomando solo una medida y después  (2)tomando varias medidas entre varios landmarks y poniéndolas todas en un un solo archivo.
(3)Finalmente les mostrare como se calculas varias distancias en varios individuos y ponerlas en un archivo como una matriz de distancias.

Primero necesitamos los archivos de los ladmarks en formato .txt, /aunque podria ser tambien .csv), llamaremos 7 archivos .txt de 7 individuos de Diachlorus curvipes a los cuales se les ha tomado los mismos landmarks (12), los archivos seran los siguientes y estaran sobre el directorio sobre el cual estemos trabajando en R.


dcp1.txt




"x" "y"
"lan1" 0.310384615384615 3.01932692307692
"lan2" 1.87769230769231 2.72125
"lan3" 1.96423076923077 3.01932692307692
"lan4" 5.06038461538461 2.58663461538461
"lan5" 5.13730769230769 2.41355769230769
"lan6" 4.40653846153846 1.75009615384615
"lan7" 4.07961538461538 1.57701923076923
"lan8" 3.06038461538462 2.05778846153846
"lan9" 3.16615384615385 1.33663461538461
"lan10" 2.07961538461538 1.32701923076923
"lan11" 2.06038461538462 1.37509615384615
"lan12" 1.74307692307692 2.25009615384615


dcp2.txt




"x" "y"
"lan1" 0.719286944871943 3.44351836802884
"lan2" 2.58291646916267 3.09794468140540
"lan3" 2.74336139509498 3.48054412016707
"lan4" 6.23612401346766 2.99920934237012
"lan5" 6.2854916829853 2.78939674692018
"lan6" 5.42155746642668 2.02419786939683
"lan7" 5.08832569718265 1.85141102608511
"lan8" 3.87881779400059 2.45616497767614
"lan9" 4.06394655469172 1.57988884373812
"lan10" 2.65696797343912 1.53052117422048
"lan11" 2.64462605605971 1.60457267849693
"lan12" 2.37310387371272 2.56724223409082


dcp3.txt




"x" "y"
"lan1" 1.10064827940554 3.7794857454421
"lan2" 2.84403223991721 3.48688983598560
"lan3" 2.92937271350869 3.86482621903358
"lan4" 6.42833213075936 3.30401739257528
"lan5" 6.44052362698672 3.18210243030174
"lan6" 5.57492739484456 2.43842116043313
"lan7" 5.29452298161541 2.29212320570488
"lan8" 4.18509682492618 2.86512352839053
"lan9" 4.31920328342707 1.97514430379367
"lan10" 3.07567066823694 1.91418682265689
"lan11" 3.03909617955487 1.97514430379367
"lan12" 2.6855427889616 2.99922998689142


dcp4.txt




"x" "y"
"lan1" 0.941808934714676 3.94869970156446
"lan2" 2.78204119586783 3.50996883135577
"lan3" 2.90391088203692 3.92432576433065
"lan4" 6.37719693785579 3.24185552178378
"lan5" 6.41375784370651 3.05905099253016
"lan6" 5.45098732297076 2.27908500104803
"lan7" 5.23162188786641 2.18158925211276
"lan8" 4.0494859320263 2.76656374572436
"lan9" 4.18354258681229 1.91347594254078
"lan10" 2.86734997618619 1.98659775424223
"lan11" 2.84297603895237 2.04753259732677
"lan12" 2.50174091767894 3.03467705529634

dcp5.txt

"x" "y"
"lan1" 0.585132420157185 3.30044582245595
"lan2" 2.00234100248879 2.90246259043132
"lan3" 2.09941008347041 3.18396292527801
"lan4" 4.87558579954464 2.49477245030853
"lan5" 4.90470652383913 2.34916882883611
"lan6" 4.19610223267333 1.79587506724089
"lan7" 3.93401571402296 1.68909907816111
"lan8" 3.02156635279576 2.27151356405081
"lan9" 3.10892852567922 1.56290927288501
"lan10" 2.03146172678327 1.65997835386663
"lan11" 2.02175481868511 1.68909907816111
"lan12" 1.78878902432923 2.53360008270118

dcp6.txt

"x" "y"
"lan1" 0.508475803224269 4.05301396510432
"lan2" 2.30040803313004 3.47574720647696
"lan3" 2.45675111359162 3.80045975820485
"lan4" 5.94440444696526 2.65795263175486
"lan5" 5.93237805616053 2.36931925244118
"lan6" 5.07850430902422 1.8521844478375
"lan7" 4.72973897568686 1.75597332139961
"lan8" 3.64736380326055 2.56174150531697
"lan9" 3.67141658487002 1.68381497657119
"lan10" 2.30040803313004 1.90029001105645
"lan11" 2.27635525152057 1.99650113749434
"lan12" 2.07190660784004 2.98266518348275

dcp7.txt

"x" "y"
"lan1" 0.137282309611791 3.15335927001144
"lan2" 1.72352748369761 2.92263269923532
"lan3" 1.87734519754835 3.2687225553995
"lan4" 5.06906275995132 2.84572384230995
"lan5" 5.12674440264535 2.68229252134353
"lan6" 4.32881501204461 1.99011280901518
"lan7" 4.07886122703715 1.86513591651145
"lan8" 3.0117508371976 2.40349791498906
"lan9" 3.18479576527969 1.62479573861966
"lan10" 1.95425405447372 1.51904606034727
"lan11" 1.92541323312671 1.59595491727264
"lan12" 1.61777780542522 2.45156595056741


Entonces...

>dcp1<-read.table ("dcp1.txt")

#llamo el archivo "dcp1.txt" de mi directorio y lo leo en R, poniéndole como nombre "dcp1" (Diachlorus curvipes, individuo 1).

> dcp1

#verifico que contiene mi nuevo objeto que acabo de llamar.

              x        y
lan1  0.3103846 3.019327
lan2  1.8776923 2.721250
lan3  1.9642308 3.019327
lan4  5.0603846 2.586635
lan5  5.1373077 2.413558
lan6  4.4065385 1.750096
lan7  4.0796154 1.577019
lan8  3.0603846 2.057788
lan9  3.1661538 1.336635
lan10 2.0796154 1.327019
lan11 2.0603846 1.375096
lan12 1.7430769 2.250096

>dist1_5<-list((sqrt(sum((dcp1[1,]-dcp1[5,])^2))))

#tomo la distancia entre el landmark 1  ("dcp1[1,]") de coordenadas  (x)0.3103846 (y)3.019327 y 5  ("dcp1[5,]") de coordenadas (X)5.1373077 (y)2.413558.

>dist1_5

#miro el valor de distancia que obtuve entre los dos puntos:

[[1]]
[1] 4.864786

> disdcp1<-list((sqrt(sum((dcp1[1,]-dcp1[5,])^2))),(sqrt(sum((dcp1[2,]-dcp1[8,])^2))),(sqrt(sum((dcp1[1,]-dcp1[11,])^2))),(sqrt(sum((dcp1[6,]-dcp1[7,])^2))),(sqrt(sum((dcp1[3,]-dcp1[10,])^2))),(sqrt(sum((dcp1[2,]-dcp1[12,])^2))),(sqrt(sum((dcp1[3,]-dcp1[4,])^2))),(sqrt(sum((dcp1[5,]-dcp1[10,])^2))),(sqrt(sum((dcp1[4,]-dcp1[6,])^2))))

#ahora tomo 9 distancias a partir de 12 landmarks, tomados de una foto. las distancias serian entre los landmarks 1y5, 2y8, 1y11, 6y7, 3y10, 2y12, 3y4, 5y10,   4y6. Y nombor este objeto "disdcp1".

> disdcp1

#reviso las distancias que obtuve:

[[1]]
[1] 4.864786

[[2]]
[1] 1.356076

[[3]]
[1] 2.401249

[[4]]
[1] 0.3699112

[[5]]
[1] 1.696237

[[6]]
[1] 0.4900074

[[7]]
[1] 3.126242

[[8]]
[1] 3.245004

[[9]]
[1] 1.061749

> distdcp1<-matrix(disdcp1, 1, 9)

#convierto mis resultados de la forma de lista a la forma de matriz!!!
> distdcp1

#reviso como se organizaron mis resultados en la nueva forma que los organice (matriz):

     [,1]     [,2]     [,3]     [,4]      [,5]     [,6]      [,7]     [,8]    
[1,] 4.864786 1.356076 2.401249 0.3699112 1.696237 0.4900074 3.126242 3.245004
     [,9]    
[1,] 1.061749

>write.table(distdcp1, file = "distdcp1.txt")

#escribo un archivo en formato .txt que contenga mi matriz de distancias halladas, que quedara sobre el directorio el cual estemos trabajando, al abrir este archivo en un edito de texto (geany, wordpad, blog de notas, etc) tendrá esta apariencia:

"V1" "V2" "V3" "V4" "V5" "V6" "V7" "V8" "V9"
"1" 4.86478599236177 1.35607607002421 2.40124859656080 0.36991123195522 1.69623669778461 0.490007396393886 3.1262423565983 3.24500355576432 1.06174921263570

> dcp1.txt<-read.table ("dcp1.txt")
> dcp2.txt<-read.table ("dcp2.txt")
> dcp3.txt<-read.table ("dcp3.txt")
> dcp4.txt<-read.table ("dcp4.txt")
> dcp5.txt<-read.table ("dcp5.txt")
> dcp6.txt<-read.table ("dcp6.txt")
> dcp7.txt<-read.table ("dcp7.txt")

#leo de mi directorio 7 archivos que contienen los landmarks de 7 individuos de la misma especie (no necesariamente tienen que ser de la misma especie, valga la aclaración) 

>disdcp1al7<-list((sqrt(sum((dcp1.txt[1,]-dcp1.txt[5,])^2))),(sqrt(sum((dcp1.txt[2,]-dcp1.txt[8,])^2))),(sqrt(sum((dcp1.txt[1,]-dcp1.txt[11,])^2))),(sqrt(sum((dcp1.txt[6,]-dcp1.txt[7,])^2))),(sqrt(sum((dcp1.txt[3,]-dcp1.txt[10,])^2))),(sqrt(sum((dcp1.txt[2,]-dcp1.txt[12,])^2))),(sqrt(sum((dcp1.txt[3,]-dcp1.txt[4,])^2))),(sqrt(sum((dcp1.txt[5,]-dcp1.txt[10,])^2))),(sqrt(sum((dcp1.txt[4,]-dcp1.txt[6,])^2))),(sqrt(sum((dcp2.txt[1,]-dcp2.txt[5,])^2))),(sqrt(sum((dcp2.txt[2,]-dcp2.txt[8,])^2))),(sqrt(sum((dcp2.txt[1,]-dcp2.txt[11,])^2))),(sqrt(sum((dcp2.txt[6,]-dcp2.txt[7,])^2))),(sqrt(sum((dcp2.txt[3,]-dcp2.txt[10,])^2))),(sqrt(sum((dcp2.txt[2,]-dcp2.txt[12,])^2))),(sqrt(sum((dcp2.txt[3,]-dcp2.txt[4,])^2))),(sqrt(sum((dcp2.txt[5,]-dcp2.txt[10,])^2))),(sqrt(sum((dcp2.txt[4,]-dcp2.txt[6,])^2))),(sqrt(sum((dcp3.txt[1,]-dcp3.txt[5,])^2))),(sqrt(sum((dcp3.txt[2,]-dcp3.txt[8,])^2))),(sqrt(sum((dcp3.txt[1,]-dcp3.txt[11,])^2))),(sqrt(sum((dcp3.txt[6,]-dcp3.txt[7,])^2))),(sqrt(sum((dcp3.txt[3,]-dcp3.txt[10,])^2))),(sqrt(sum((dcp3.txt[2,]-dcp3.txt[12,])^2))),(sqrt(sum((dcp3.txt[3,]-dcp3.txt[4,])^2))),(sqrt(sum((dcp3.txt[5,]-dcp3.txt[10,])^2))),(sqrt(sum((dcp3.txt[4,]-dcp3.txt[6,])^2))),(sqrt(sum((dcp4.txt[1,]-dcp4.txt[5,])^2))),(sqrt(sum((dcp4.txt[2,]-dcp4.txt[8,])^2))),(sqrt(sum((dcp4.txt[1,]-dcp4.txt[11,])^2))),(sqrt(sum((dcp4.txt[6,]-dcp4.txt[7,])^2))),(sqrt(sum((dcp4.txt[3,]-dcp4.txt[10,])^2))),(sqrt(sum((dcp4.txt[2,]-dcp4.txt[12,])^2))),(sqrt(sum((dcp4.txt[3,]-dcp4.txt[4,])^2))),(sqrt(sum((dcp4.txt[5,]-dcp4.txt[10,])^2))),(sqrt(sum((dcp4.txt[4,]-dcp4.txt[6,])^2))),(sqrt(sum((dcp5.txt[1,]-dcp5.txt[5,])^2))),(sqrt(sum((dcp5.txt[2,]-dcp5.txt[8,])^2))),(sqrt(sum((dcp5.txt[1,]-dcp5.txt[11,])^2))),(sqrt(sum((dcp5.txt[6,]-dcp5.txt[7,])^2))),(sqrt(sum((dcp5.txt[3,]-dcp5.txt[10,])^2))),(sqrt(sum((dcp5.txt[2,]-dcp5.txt[12,])^2))),(sqrt(sum((dcp5.txt[3,]-dcp5.txt[4,])^2))),(sqrt(sum((dcp5.txt[5,]-dcp5.txt[10,])^2))),(sqrt(sum((dcp5.txt[4,]-dcp5.txt[6,])^2))),(sqrt(sum((dcp6.txt[1,]-dcp6.txt[5,])^2))),(sqrt(sum((dcp6.txt[2,]-dcp6.txt[8,])^2))),(sqrt(sum((dcp6.txt[1,]-dcp6.txt[11,])^2))),(sqrt(sum((dcp6.txt[6,]-dcp6.txt[7,])^2))),(sqrt(sum((dcp6.txt[3,]-dcp6.txt[10,])^2))),(sqrt(sum((dcp6.txt[2,]-dcp6.txt[12,])^2))),(sqrt(sum((dcp6.txt[3,]-dcp6.txt[4,])^2))),(sqrt(sum((dcp6.txt[5,]-dcp6.txt[10,])^2))),(sqrt(sum((dcp6.txt[4,]-dcp6.txt[6,])^2))),(sqrt(sum((dcp7.txt[1,]-dcp7.txt[5,])^2))),(sqrt(sum((dcp7.txt[2,]-dcp7.txt[8,])^2))),(sqrt(sum((dcp7.txt[1,]-dcp7.txt[11,])^2))),(sqrt(sum((dcp7.txt[6,]-dcp7.txt[7,])^2))),(sqrt(sum((dcp7.txt[3,]-dcp7.txt[10,])^2))),(sqrt(sum((dcp7.txt[2,]-dcp7.txt[12,])^2))),(sqrt(sum((dcp7.txt[3,]-dcp7.txt[4,])^2))),(sqrt(sum((dcp7.txt[5,]-dcp7.txt[10,])^2))),(sqrt(sum((dcp7.txt[4,]-dcp7.txt[6,])^2))))

#Tomo las mismas 9 medidas de cada uno de los siete individuos que tengo de la forma como lo hacíamos, con un solo individuo las nueve medidas, pero esta vez se hacen con los 7 individuos poniendo los comandos uno tras de otro. 


>disdcp1al7

#reviso las medidas que obtuve:

[[1]]
[1] 4.864786

[[2]]
[1] 1.356076

[[3]]
[1] 2.401249

[[4]]
[1] 0.3699112

[[5]]
[1] 1.696237

[[6]]
[1] 0.4900074

[[7]]
[1] 3.126242

[[8]]
[1] 3.245004

[[9]]
[1] 1.061749

[[10]]
[1] 5.604508

[[11]]
[1] 1.446113

[[12]]
[1] 2.662452

[[13]]
[1] 0.3753648

[[14]]
[1] 1.951936

[[15]]
[1] 0.5706719

[[16]]
[1] 3.525773

[[17]]
[1] 3.840697

[[18]]
[1] 1.270498

[[19]]
[1] 5.373187

[[20]]
[1] 1.478191

[[21]]
[1] 2.64825

[[22]]
[1] 0.3162748

[[23]]
[1] 1.956118

[[24]]
[1] 0.512768

[[25]]
[1] 3.543617

[[26]]
[1] 3.595809

[[27]]
[1] 1.215548

[[28]]
[1] 5.543798

[[29]]
[1] 1.469376

[[30]]
[1] 2.688656

[[31]]
[1] 0.2400554

[[32]]
[1] 1.938073

[[33]]
[1] 0.5517885

[[34]]
[1] 3.539701

[[35]]
[1] 3.705019

[[36]]
[1] 1.335961

[[37]]
[1] 4.423081

[[38]]
[1] 1.198715

[[39]]
[1] 2.158778

[[40]]
[1] 0.2830026

[[41]]
[1] 1.525499

[[42]]
[1] 0.4262206

[[43]]
[1] 2.860443

[[44]]
[1] 2.954745

[[45]]
[1] 0.9747592

[[46]]
[1] 5.67922

[[47]]
[1] 1.627789

[[48]]
[1] 2.711944

[[49]]
[1] 0.3617925

[[50]]
[1] 1.906591

[[51]]
[1] 0.5434545

[[52]]
[1] 3.67002

[[53]]
[1] 3.66213

[[54]]
[1] 1.182813

[[55]]
[1] 5.01165

[[56]]
[1] 1.388892

[[57]]
[1] 2.37127

[[58]]
[1] 0.2794568

[[59]]
[1] 1.751366

[[60]]
[1] 0.4827907

[[61]]
[1] 3.219626

[[62]]
[1] 3.379029

[[63]]
[1] 1.131387 


> distdcp1al7<-matrix(disdcp1al7, 9, 7)

#organizo mis medidas en una matriz de 9x7 (9 medidas correspondientes a 7 individuos).

> distdcp1al7

#ahora mis resultados se verán asi:

[,1]      [,2]      [,3]      [,4]      [,5]      [,6]      [,7]     
 [1,] 4.864786  5.604508  5.373187  5.543798  4.423081  5.67922   5.01165  
 [2,] 1.356076  1.446113  1.478191  1.469376  1.198715  1.627789  1.388892 
 [3,] 2.401249  2.662452  2.64825   2.688656  2.158778  2.711944  2.37127  
 [4,] 0.3699112 0.3753648 0.3162748 0.2400554 0.2830026 0.3617925 0.2794568
 [5,] 1.696237  1.951936  1.956118  1.938073  1.525499  1.906591  1.751366 
 [6,] 0.4900074 0.5706719 0.512768  0.5517885 0.4262206 0.5434545 0.4827907
 [7,] 3.126242  3.525773  3.543617  3.539701  2.860443  3.67002   3.219626 
 [8,] 3.245004  3.840697  3.595809  3.705019  2.954745  3.66213   3.379029 
 [9,] 1.061749  1.270498  1.215548  1.335961  0.9747592 1.182813  1.131387 


>write.table(distdcp1al7, file = "distdcp1al7.txt")

#por ultimo escribo mis resultados en una tabla en formato .txt para depues realizar los correspondientes analisis que desee. Este archivo se vera asi:

"V1" "V2" "V3" "V4" "V5" "V6" "V7"
"1" 4.86478599236177 5.60450803209322 5.37318672231976 5.54379830875036 4.42308131914733 5.67922032807822 5.01165001367009
"2" 1.35607607002421 1.44611252391904 1.47819063805238 1.46937649380901 1.19871472360978 1.62778876652987 1.38889176421965
"3" 2.40124859656080 2.66245224222281 2.64825004477358 2.68865630315053 2.15877795206031 2.71194445994619 2.37126981121088
"4" 0.36991123195522 0.375364762937514 0.316274764273176 0.240055441885765 0.283002570840330 0.361792535285579 0.279456827255025
"5" 1.69623669778461 1.95193578609124 1.95611787637005 1.93807289360278 1.52549859163370 1.90659073395059 1.75136598391422
"6" 0.490007396393886 0.570671895921929 0.512768012343564 0.551788472462498 0.426220596622018 0.54345449005379 0.48279071672441
"7" 3.1262423565983 3.52577280558723 3.54361729640782 3.53970078643772 2.86044316799144 3.67002020536417 3.21962558528581
"8" 3.24500355576432 3.84069681415740 3.59580942516654 3.70501885435476 2.95474519622204 3.66212980073973 3.37902905259225
"9" 1.06174921263570 1.27049834001773 1.21554781081638 1.33596082513795 0.97475918549504 1.18281250208974 1.13138718773756



De este modo, se pueden llegar a tomar muchísimas distancias de muchísimos individuos, haciendo los scripts y corriéndolos en R, y se obtienen mucho mas fácil,  confiables y rápido. 

que les rinda!!!!



No hay comentarios:

Publicar un comentario