Regressão linear bayesiana - Bayesian linear regression

Em estatística , a regressão linear bayesiana é uma abordagem da regressão linear em que a análise estatística é realizada no contexto da inferência bayesiana . Quando o modelo de regressão tem erros com distribuição normal , e se uma forma particular de distribuição anterior é assumida, resultados explícitos estão disponíveis para as distribuições de probabilidade posteriores dos parâmetros do modelo.

Configuração do modelo

Considere um problema de regressão linear padrão , em que especificamos a média da distribuição condicional de um determinado vetor preditor : ${\ displaystyle i = 1, \ ldots, n}$ ${\ displaystyle y_ {i}}$ ${\ displaystyle k \ times 1}$ ${\ displaystyle \ mathbf {x} _ {i}}$

{\ displaystyle y_ {i} = \ mathbf {x} _ {i} ^ {\ rm {T}} {\ boldsymbol {\ beta}} + \ varepsilon _ {i},}

onde é um vetor, e são variáveis aleatórias independentes e normalmente distribuídas de forma idêntica : ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle k \ times 1}$ ${\ displaystyle \ varepsilon _ {i}}$

{\ displaystyle \ varepsilon _ {i} \ sim N (0, \ sigma ^ {2}).}

Isso corresponde à seguinte função de verossimilhança :

{\ displaystyle \ rho (\ mathbf {y} \ mid \ mathbf {X}, {\ boldsymbol {\ beta}}, \ sigma ^ {2}) \ propto (\ sigma ^ {2}) ^ {- n / 2} \ exp \ left (- {\ frac {1} {2 \ sigma ^ {2}}} (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta}}) ^ {\ rm { T}} (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta}}) \ right).}

A solução de mínimos quadrados ordinários é usada para estimar o vetor de coeficiente usando o pseudoinverso Moore-Penrose :

{\ displaystyle {\ hat {\ boldsymbol {\ beta}}} = (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X}) ^ {- 1} \ mathbf {X} ^ {\ rm {T}} \ mathbf {y}}

onde está a matriz de design , cada linha da qual é um vetor preditor ; e é a coluna -vetor . ${\ displaystyle \ mathbf {X}}$ ${\ displaystyle n \ times (k + 1)}$ ${\ displaystyle \ mathbf {x} _ {i} ^ {\ rm {T}}}$ ${\ displaystyle \ mathbf {y}}$ ${\ displaystyle n}$ ${\ displaystyle [y_ {1} \; \ cdots \; y_ {n}] ^ {\ rm {T}}}$

Esta é uma abordagem frequentista e assume que há medidas suficientes para dizer algo significativo . Na abordagem bayesiana , os dados são complementados com informações adicionais na forma de uma distribuição de probabilidade anterior . A crença anterior sobre os parâmetros é combinada com a função de verossimilhança dos dados de acordo com o teorema de Bayes para produzir a crença posterior sobre os parâmetros e . O prior pode assumir diferentes formas funcionais dependendo do domínio e das informações disponíveis a priori . ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle \ sigma}$

Com antecedentes conjugados

Conjugado distribuição anterior

Para uma distribuição anterior arbitrária, pode não haver solução analítica para a distribuição posterior . Nesta seção, consideraremos um assim chamado conjugado anterior para o qual a distribuição posterior pode ser derivada analiticamente.

Um prior é conjugado a esta função de verossimilhança se tiver a mesma forma funcional em relação a e . Uma vez que a probabilidade logarítmica é quadrática em , a probabilidade logarítmica é reescrita de modo que a probabilidade se torne normal em . Escreva ${\ displaystyle \ rho ({\ boldsymbol {\ beta}}, \ sigma ^ {2})}$ ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle \ sigma}$ ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle ({\ boldsymbol {\ beta}} - {\ hat {\ boldsymbol {\ beta}}})}$

{\ displaystyle (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta}}) ^ {\ rm {T}} (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta }}) = (\ mathbf {y} - \ mathbf {X} {\ hat {\ boldsymbol {\ beta}}}) ^ {\ rm {T}} (\ mathbf {y} - \ mathbf {X} { \ hat {\ boldsymbol {\ beta}}}) + ({\ boldsymbol {\ beta}} - {\ hat {\ boldsymbol {\ beta}}}) ^ {\ rm {T}} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X}) ({\ boldsymbol {\ beta}} - {\ hat {\ boldsymbol {\ beta}}}).}

A probabilidade agora foi reescrita como

{\ displaystyle \ rho (\ mathbf {y} | \ mathbf {X}, {\ boldsymbol {\ beta}}, \ sigma ^ {2}) \ propto (\ sigma ^ {2}) ^ {- {\ frac {v} {2}}} \ exp \ left (- {\ frac {vs ^ {2}} {2 {\ sigma} ^ {2}}} \ right) (\ sigma ^ {2}) ^ {- {\ frac {nv} {2}}} \ exp \ left (- {\ frac {1} {2 {\ sigma} ^ {2}}} ({\ boldsymbol {\ beta}} - {\ hat {\ boldsymbol {\ beta}}}) ^ {\ rm {T}} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X}) ({\ boldsymbol {\ beta}} - {\ hat { \ boldsymbol {\ beta}}}) \ right),}

Onde

{\ displaystyle vs ^ {2} = (\ mathbf {y} - \ mathbf {X} {\ hat {\ boldsymbol {\ beta}}}) ^ {\ rm {T}} (\ mathbf {y} - \ mathbf {X} {\ hat {\ boldsymbol {\ beta}}}) \ quad {\ text {e}} \ quad v = nk,}

onde é o número de coeficientes de regressão. ${\ displaystyle k}$

Isso sugere um formulário para o anterior:

{\ displaystyle \ rho ({\ boldsymbol {\ beta}}, \ sigma ^ {2}) = \ rho (\ sigma ^ {2}) \ rho ({\ boldsymbol {\ beta}} \ mid \ sigma ^ { 2}),}

onde é uma distribuição gama inversa ${\ displaystyle \ rho (\ sigma ^ {2})}$

{\ displaystyle \ rho (\ sigma ^ {2}) \ propto (\ sigma ^ {2}) ^ {- {\ frac {v_ {0}} {2}} - 1} \ exp \ left (- {\ frac {v_ {0} s_ {0} ^ {2}} {2 \ sigma ^ {2}}} \ right).}

Na notação introduzida no artigo de distribuição gama inversa , esta é a densidade de uma distribuição com e com e como os valores anteriores de e , respectivamente. Equivalentemente, também pode ser descrito como uma distribuição qui-quadrada inversa em escala , ${\ displaystyle {\ text {Inv-Gamma}} (a_ {0}, b_ {0})}$ ${\ displaystyle a_ {0} = {\ tfrac {v_ {0}} {2}}}$ ${\ displaystyle b_ {0} = {\ tfrac {1} {2}} v_ {0} s_ {0} ^ {2}}$ ${\ displaystyle v_ {0}}$ ${\ displaystyle s_ {0} ^ {2}}$ ${\ displaystyle v}$ ${\ displaystyle s ^ {2}}$ ${\ displaystyle {\ text {escala-inv -}} \ chi ^ {2} (v_ {0}, s_ {0} ^ {2}).}$

Além disso, a densidade prévia condicional é uma distribuição normal , ${\ displaystyle \ rho ({\ boldsymbol {\ beta}} | \ sigma ^ {2})}$

{\ displaystyle \ rho ({\ boldsymbol {\ beta}} \ mid \ sigma ^ {2}) \ propto (\ sigma ^ {2}) ^ {- k / 2} \ exp \ left (- {\ frac { 1} {2 \ sigma ^ {2}}} ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {0}) ^ {\ rm {T}} \ mathbf {\ Lambda} _ {0} ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {0}) \ right).}

Na notação da distribuição normal , a distribuição condicional a priori é ${\ displaystyle {\ mathcal {N}} \ left ({\ boldsymbol {\ mu}} _ {0}, \ sigma ^ {2} \ mathbf {\ Lambda} _ {0} ^ {- 1} \ right) .}$

Distribuição posterior

Com o prior agora especificado, a distribuição posterior pode ser expressa como

{\ displaystyle {\ begin {alinhados} \ rho ({\ boldsymbol {\ beta}}, \ sigma ^ {2} \ mid \ mathbf {y}, \ mathbf {X}) & \ propto \ rho (\ mathbf { y} \ mid \ mathbf {X}, {\ boldsymbol {\ beta}}, \ sigma ^ {2}) \ rho ({\ boldsymbol {\ beta}} \ mid \ sigma ^ {2}) \ rho (\ sigma ^ {2}) \\ & \ propto (\ sigma ^ {2}) ^ {- n / 2} \ exp \ left (- {\ frac {1} {2 {\ sigma} ^ {2}}} (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta}}) ^ {\ rm {T}} (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta}}) \ right) (\ sigma ^ {2}) ^ {- k / 2} \ exp \ left (- {\ frac {1} {2 \ sigma ^ {2}}} ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {0}) ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {0} ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu} } _ {0}) \ right) (\ sigma ^ {2}) ^ {- (a_ {0} +1)} \ exp \ left (- {\ frac {b_ {0}} {\ sigma ^ {2 }}} \ right) \ end {alinhados}}}

Com algum rearranjo, o posterior pode ser reescrito de modo que a média posterior do vetor de parâmetro possa ser expressa em termos do estimador de mínimos quadrados e a média anterior , com a força do anterior indicada pela matriz de precisão anterior ${\ displaystyle {\ boldsymbol {\ mu}} _ {n}}$ ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle {\ hat {\ boldsymbol {\ beta}}}}$ ${\ displaystyle {\ boldsymbol {\ mu}} _ {0}}$ ${\ displaystyle {\ boldsymbol {\ Lambda}} _ {0}}$

{\ displaystyle {\ boldsymbol {\ mu}} _ {n} = (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + {\ boldsymbol {\ Lambda}} _ {0}) ^ {-1} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} {\ hat {\ boldsymbol {\ beta}}} + {\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0}).}

Para justificar que é de fato a média posterior, os termos quadráticos no exponencial podem ser reorganizados como uma forma quadrática em . ${\ displaystyle {\ boldsymbol {\ mu}} _ {n}}$ ${\ displaystyle {\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {n}}$

{\ displaystyle (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta}}) ^ {\ rm {T}} (\ mathbf {y} - \ mathbf {X} {\ boldsymbol {\ beta }}) + ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {0}) ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {0} ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {0}) = ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {n}) ^ {\ rm {T }} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + {\ boldsymbol {\ Lambda}} _ {0}) ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {n}) + \ mathbf {y} ^ {\ rm {T}} \ mathbf {y} - {\ boldsymbol {\ mu}} _ {n} ^ {\ rm {T}} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + {\ boldsymbol {\ Lambda}} _ {0}) {\ boldsymbol {\ mu}} _ {n} + {\ boldsymbol {\ mu }} _ {0} ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0}.}

Agora, o posterior pode ser expresso como uma distribuição normal vezes uma distribuição gama inversa :

{\ displaystyle \ rho ({\ boldsymbol {\ beta}}, \ sigma ^ {2} \ mid \ mathbf {y}, \ mathbf {X}) \ propto (\ sigma ^ {2}) ^ {- k / 2} \ exp \ left (- {\ frac {1} {2 {\ sigma} ^ {2}}} ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {n}) ^ {\ rm {T}} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + \ mathbf {\ Lambda} _ {0}) ({\ boldsymbol {\ beta}} - {\ boldsymbol {\ mu}} _ {n}) \ right) (\ sigma ^ {2}) ^ {- {\ frac {n + 2a_ {0}} {2}} - 1} \ exp \ left (- { \ frac {2b_ {0} + \ mathbf {y} ^ {\ rm {T}} \ mathbf {y} - {\ boldsymbol {\ mu}} _ {n} ^ {\ rm {T}} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + {\ boldsymbol {\ Lambda}} _ {0}) {\ boldsymbol {\ mu}} _ {n} + {\ boldsymbol {\ mu} } _ {0} ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0}} {2 \ sigma ^ {2}}} \ right ).}

Portanto, a distribuição posterior pode ser parametrizada da seguinte maneira.

{\ displaystyle \ rho ({\ boldsymbol {\ beta}}, \ sigma ^ {2} \ mid \ mathbf {y}, \ mathbf {X}) \ propto \ rho ({\ boldsymbol {\ beta}} \ mid \ sigma ^ {2}, \ mathbf {y}, \ mathbf {X}) \ rho (\ sigma ^ {2} \ mid \ mathbf {y}, \ mathbf {X}),}

onde os dois fatores correspondem às densidades e distribuições, com os parâmetros destes dados por ${\ displaystyle {\ mathcal {N}} \ left ({\ boldsymbol {\ mu}} _ {n}, \ sigma ^ {2} {\ boldsymbol {\ Lambda}} _ {n} ^ {- 1} \ direito)\,}$ ${\ displaystyle {\ text {Inv-Gamma}} \ left (a_ {n}, b_ {n} \ right)}$

{\ displaystyle {\ boldsymbol {\ Lambda}} _ {n} = (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + \ mathbf {\ Lambda} _ {0}), \ quad {\ boldsymbol {\ mu}} _ {n} = ({\ boldsymbol {\ Lambda}} _ {n}) ^ {- 1} (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X } {\ hat {\ boldsymbol {\ beta}}} + {\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0}),}

{\ displaystyle a_ {n} = a_ {0} + {\ frac {n} {2}}, \ qquad b_ {n} = b_ {0} + {\ frac {1} {2}} (\ mathbf { y} ^ {\ rm {T}} \ mathbf {y} + {\ boldsymbol {\ mu}} _ {0} ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0} - {\ boldsymbol {\ mu}} _ {n} ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {n} {\ boldsymbol {\ mu} } _ {n}).}

Isso pode ser interpretado como aprendizagem Bayesiana, onde os parâmetros são atualizados de acordo com as seguintes equações.

{\ displaystyle {\ boldsymbol {\ mu}} _ {n} = (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + {\ boldsymbol {\ Lambda}} _ {0}) ^ {-1} ({\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0} + \ mathbf {X} ^ {\ rm {T}} \ mathbf {X} {\ hat {\ boldsymbol {\ beta}}}),}

{\ displaystyle {\ boldsymbol {\ Lambda}} _ {n} = (\ mathbf {X} ^ {\ rm {T}} \ mathbf {X} + {\ boldsymbol {\ Lambda}} _ {0}), }

{\ displaystyle a_ {n} = a_ {0} + {\ frac {n} {2}},}

{\ displaystyle b_ {n} = b_ {0} + {\ frac {1} {2}} (\ mathbf {y} ^ {\ rm {T}} \ mathbf {y} + {\ boldsymbol {\ mu} } _ {0} ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {0} {\ boldsymbol {\ mu}} _ {0} - {\ boldsymbol {\ mu}} _ {n} ^ {\ rm {T}} {\ boldsymbol {\ Lambda}} _ {n} {\ boldsymbol {\ mu}} _ {n}).}

Evidência de modelo

A evidência do modelo é a probabilidade dos dados dados o modelo . Também é conhecido como probabilidade marginal e como densidade preditiva anterior . Aqui, o modelo é definido pela função de verossimilhança e pela distribuição anterior dos parâmetros, ou seja . A evidência do modelo captura em um único número quão bem esse modelo explica as observações. A evidência do modelo do modelo de regressão linear Bayesiana apresentada nesta seção pode ser usada para comparar modelos lineares concorrentes por comparação de modelo Bayesiano . Esses modelos podem diferir no número e nos valores das variáveis preditoras, bem como em seus antecedentes nos parâmetros do modelo. A complexidade do modelo já é levada em consideração pela evidência do modelo, porque ela marginaliza os parâmetros ao integrar todos os valores possíveis de e . ${\ displaystyle p (\ mathbf {y} \ mid m)}$ ${\ displaystyle m}$ ${\ displaystyle p (\ mathbf {y} \ mid \ mathbf {X}, {\ boldsymbol {\ beta}}, \ sigma)}$ ${\ displaystyle p ({\ boldsymbol {\ beta}}, \ sigma)}$ ${\ displaystyle p (\ mathbf {y}, {\ boldsymbol {\ beta}}, \ sigma \ mid \ mathbf {X})}$ ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle \ sigma}$

{\ displaystyle p (\ mathbf {y} | m) = \ int p (\ mathbf {y} \ mid \ mathbf {X}, {\ boldsymbol {\ beta}}, \ sigma) \, p ({\ boldsymbol {\ beta}}, \ sigma) \, d {\ boldsymbol {\ beta}} \, d \ sigma}

Essa integral pode ser calculada analiticamente e a solução é dada na seguinte equação.

{\ displaystyle p (\ mathbf {y} \ mid m) = {\ frac {1} {(2 \ pi) ^ {n / 2}}} {\ sqrt {\ frac {\ det ({\ boldsymbol {\ Lambda}} _ {0})} {\ det ({\ boldsymbol {\ Lambda}} _ {n})}}} \ cdot {\ frac {b_ {0} ^ {a_ {0}}} {b_ { n} ^ {a_ {n}}}} \ cdot {\ frac {\ Gamma (a_ {n})} {\ Gamma (a_ {0})}}}

Aqui denota a função gama . Como escolhemos um conjugado a priori, a probabilidade marginal também pode ser facilmente calculada avaliando a seguinte igualdade para valores arbitrários de e . ${\ displaystyle \ Gamma}$ ${\ displaystyle {\ boldsymbol {\ beta}}}$ ${\ displaystyle \ sigma}$

{\ displaystyle p (\ mathbf {y} \ mid m) = {\ frac {p ({\ boldsymbol {\ beta}}, \ sigma | m) \, p (\ mathbf {y} \ mid \ mathbf {X }, {\ boldsymbol {\ beta}}, \ sigma, m)} {p ({\ boldsymbol {\ beta}}, \ sigma \ mid \ mathbf {y}, \ mathbf {X}, m)}}}

Observe que esta equação nada mais é do que um rearranjo do teorema de Bayes . Inserir as fórmulas para o anterior, a verossimilhança e o posterior e simplificar a expressão resultante leva à expressão analítica fornecida acima.

Outros casos

Em geral, pode ser impossível ou impraticável derivar a distribuição posterior analiticamente. No entanto, é possível aproximar a posterior por um método de inferência Bayesiana aproximada , como a amostragem de Monte Carlo ou Bayes variacional .

O caso especial é chamado de regressão de crista . ${\ displaystyle {\ boldsymbol {\ mu}} _ {0} = 0, \ mathbf {\ Lambda} _ {0} = c \ mathbf {I}}$

Uma análise semelhante pode ser realizada para o caso geral da regressão multivariada e parte disso fornece a estimativa bayesiana de matrizes de covariância : ver regressão linear multivariada bayesiana .

Veja também

Notas

Referências

Box, GEP ; Tiao, GC (1973). Inferência Bayesiana em Análise Estatística . Wiley. ISBN 0-471-57428-7 .
Carlin, Bradley P .; Louis, Thomas A. (2008). Métodos Bayesianos para Análise de Dados, Terceira Edição . Boca Raton, FL: Chapman e Hall / CRC. ISBN 1-58488-697-8 .
Fahrmeir, L .; Kneib, T .; Lang, S. (2009). Regressão. Modelle, Methoden und Anwendungen (Segunda ed.). Heidelberg: Springer. doi : 10.1007 / 978-3-642-01837-4 . ISBN 978-3-642-01836-7 .
Fornalski KW; Parzych G .; Pylak M .; Satuła D .; Dobrzyński L. (2010). "Aplicação do raciocínio Bayesiano e do Método da Máxima Entropia a alguns problemas de reconstrução" . Acta Physica Polonica Uma . 117 (6): 892–899. doi : 10.12693 / APhysPolA.117.892 .
Fornalski, Krzysztof W. (2015). "Aplicações da análise de regressão Bayesiana robusta". International Journal of Society Systems Science . 7 (4): 314–333. doi : 10.1504 / IJSSS.2015.073223 .
Gelman, Andrew ; Carlin, John B .; Stern, Hal S .; Rubin, Donald B. (2003). Análise de dados bayesiana, segunda edição . Boca Raton, FL: Chapman e Hall / CRC. ISBN 1-58488-388-X .
Goldstein, Michael; Wooff, David (2007). Estatística Linear de Bayes, Teoria e Métodos . Wiley. ISBN 978-0-470-01562-9 .
Minka, Thomas P. (2001) Bayesian Linear Regression , página de pesquisa da Microsoft
Rossi, Peter E .; Allenby, Greg M .; McCulloch, Robert (2006). Estatística Bayesiana e Marketing . John Wiley & Sons. ISBN 0470863676 .
O'Hagan, Anthony (1994). Inferência Bayesiana . Teoria Avançada de Estatística de Kendall. 2B (primeira edição). Halsted. ISBN 0-340-52922-9 .
Sivia, DS; Skilling, J. (2006). Análise de dados - um tutorial bayesiano (segunda edição). Imprensa da Universidade de Oxford.
Walter, Gero; Augustin, Thomas (2009). "Regressão Linear Bayesiana - Diferentes Modelos Conjugados e Sua (In) Sensibilidade ao Conflito de Dados Prévios" (PDF) . Relatório Técnico Número 069, Departamento de Estatística, Universidade de Munique .

links externos

Estimativa bayesiana de modelos lineares (wikibook de programação R) . Regressão linear Bayesian como implementado em R .

Languages

In other projects