Regressão para um modelo de forma ?

22

Eu tenho um conjunto de dados que é estatísticas de um fórum de discussão na web. Eu estou olhando para a distribuição do número de respostas que um tópico deve ter. Em particular, criei um conjunto de dados com uma lista de contagens de respostas de tópicos e, em seguida, a contagem de tópicos com esse número de respostas.

"num_replies","count"
0,627568
1,156371
2,151670
3,79094
4,59473
5,39895
6,30947
7,23329
8,18726

Se plotar o conjunto de dados em um gráfico de log-log, obtenho o que é basicamente uma linha reta:

Dados plotados em escala log-log

(Esta é uma distribuição Zipfian ). A Wikipedia me diz que linhas retas nos gráficos de log-log implicam uma função que pode ser modelada por um monômio da forma . E, de fato, observei essa função:y=axk

lines(data$num_replies, 480000 * data$num_replies ^ -1.62, col="green")

Modelo ocular

Meus globos oculares obviamente não são tão precisos quanto R. Então, como posso fazer com que o R ajuste os parâmetros deste modelo para mim com mais precisão? Tentei regressão polinomial, mas não acho que R tente ajustar o expoente como parâmetro - qual é o nome adequado para o modelo que eu quero?

Edit: Obrigado pelas respostas a todos. Como sugerido, agora ajustei um modelo linear nos logs dos dados de entrada, usando esta receita:

data <- read.csv(file="result.txt")

# Avoid taking the log of zero:
data$num_replies = data$num_replies + 1

plot(data$num_replies, data$count, log="xy", cex=0.8)

# Fit just the first 100 points in the series:
model <- lm(log(data$count[1:100]) ~ log(data$num_replies[1:100]))

points(data$num_replies, round(exp(coef(model)[1] + coef(model)[2] * log(data$num_replies))), 
       col="red")

O resultado é este, mostrando o modelo em vermelho:

Modelo equipado

Parece uma boa aproximação para meus propósitos.

Se eu usar esse modelo Zipfian (alpha = 1.703164) junto com um gerador de números aleatórios para gerar o mesmo número total de tópicos (1400930) que o conjunto de dados medido original contido (usando este código C que encontrei na web ), o resultado será semelhante. gostar:

Resultados gerados por números aleatórios

Os pontos medidos estão em preto, os gerados aleatoriamente de acordo com o modelo estão em vermelho.

Acho que isso mostra que a variação simples criada pela geração aleatória desses 1400930 pontos é uma boa explicação para a forma do gráfico original.

Se você estiver interessado em jogar com os dados brutos, eu os publiquei aqui .

thenickdude
fonte
2
Por que não pegar registros de ambas as contagens e num_replies e ajustar um modelo linear padrão a eles?
gung - Restabelece Monica
3
O que é esse enorme aumento nas contagens logo abaixo de 10.000 respostas?
Glen_b -Reinstate Monica
3
Nem as contagens nem as contagens de log apresentam variação constante (para contagens, a variação aumentará com a média, para as contagens de log geralmente diminuirá com a média). Dado que ambas as variáveis ​​são contagens e muitas contagens são muito pequenas, eu me inclinaria para um GLM binomial Poisson, quase-Poisson ou negativo, talvez com um link de log. Se você precisar usar regressão comum, pelo menos lide com a questão da variação. Outra alternativa é fazer uma transformação Anscombe ou Freeman-Tukey das contagens e ajustar um modelo de mínimos quadrados não linear.
Glen_b -Reinstate Monica
1
Esse pico interessante se deve a um "tamanho máximo de tópico" imposto por humanos em vários fóruns.
Thenickdude 23/05
2
Fudge é delicioso :) Mais prosaicamente, não há diferença entre (num_replies + 1) e (num_posts_in_topic).
thenickdude

Respostas:

22

Seu exemplo é muito bom, porque aponta claramente problemas recorrentes com esses dados.

Dois nomes comuns são função de poder e lei de poder. Na biologia e em alguns outros campos, as pessoas costumam falar de alometria, especialmente quando você está relacionando medidas de tamanho. Na física, e em alguns outros campos, as pessoas falam de leis de escala.

Eu não consideraria monomial um bom termo aqui, pois associo isso a potências inteiras. Pela mesma razão, isso não é melhor considerado como um caso especial de um polinômio.

Os problemas de adaptação de uma lei de potência à cauda de uma distribuição se transformam em problemas de adaptação de uma lei de potência à relação entre duas variáveis ​​diferentes.

A maneira mais fácil de ajustar uma lei de energia é usar logaritmos de ambas as variáveis ​​e depois ajustar uma linha reta usando regressão. Há muitas objeções a isso sempre que ambas as variáveis ​​estão sujeitas a erro, como é comum. O exemplo aqui é um caso em questão, pois ambas as variáveis ​​(e nenhuma) podem ser consideradas como resposta (variável dependente). Esse argumento leva a um método mais simétrico de ajuste.

Além disso, sempre há a questão de suposições sobre a estrutura de erros. Novamente, o exemplo aqui é um exemplo, pois os erros são claramente heterocedásticos. Isso sugere algo mais como mínimos quadrados ponderados.

Uma excelente revisão é http://www.ncbi.nlm.nih.gov/pubmed/16573844

Outro problema é que as pessoas geralmente identificam leis de energia apenas em alguns dados. As questões se tornam científicas e estatísticas, indo até a identificação de leis de poder como pensamento positivo ou um passatempo amador na moda. Grande parte da discussão surge sob os títulos fractal e comportamento sem escala, com discussões associadas que vão da física à metafísica. No seu exemplo específico, um pouco de curvatura parece evidente.

Os entusiastas das leis de poder nem sempre são compatíveis com os céticos, porque os entusiastas publicam mais do que os céticos. Sugiro que um gráfico de dispersão em escalas logarítmicas, embora um gráfico natural e excelente que seja essencial, seja acompanhado por gráficos residuais de algum tipo para verificar se há desvios da forma da função de potência.

Nick Cox
fonte
2
Obrigado, isso explica por que não consegui encontrar nada parecido com isso onde as pessoas discutiam "regressão polinomial". Atualizei minha pergunta com os resultados da adaptação desse modelo!
thenickdude
Se você está procurando uma abordagem um pouco mais rigorosa da lei de ajuste de potência e testes de significância para o modelo ajustado, provavelmente deseja este artigo: arxiv.org/abs/0706.1062 e o código que acompanha: tuvalu.santafe.edu/ ~ aaronc / powerlaws
Martin O'Leary
2
O artigo citado acima é para distribuições que são leis de poder, não para relacionamentos entre variáveis ​​que são leis de poder. O título desta pergunta se encaixa melhor na última; o exemplo dessa pergunta se encaixa melhor na primeira.
Nick Cox
1

Se você presumir que uma potência é um bom modelo para se ajustar, poderá usar log(y) ~ log(x)como modelo e ajustar uma regressão linear usando lm():

Tente o seguinte:

# Generate some data
set.seed(42)

x <- seq(1, 10, 1)

a = 10
b = 2
scatt <- rnorm(10, sd = 0.2)


dat <- data.frame(
  x = x,
  y = a*x^(-b) + scatt
)

Ajustar um modelo:

# Fit a model
model <- lm(log(y) ~ log(x) + 1, data = dat) 
summary(model)

pred <- data.frame(
  x = dat$x,
  p = exp(predict(model, dat))
)

Agora crie um gráfico:

# Create a plot
library(ggplot2)
ggplot() +
  geom_point(data = dat, aes(x=x, y=y)) +
  geom_line(data = pred, aes(x=x, y=p), col = "red")

insira a descrição da imagem aqui

Andrie
fonte