[julia-users] Re: [help]: configuring Julia and Sublime text 2 on windows 8.1

2015-11-04 Thread Donald Lacombe
John,

I have not tried Sublime Text with Julia 0.4 (32 or 64-bit). I was able to 
install and run an earlier version of Julia (I believe some 0.3 variety) 
64-bit and it worked just fine.

If you look at the issues related to Sublime-IJulia it appears that there 
are some issues with 0.4: https://github.com/quinnj/Sublime-IJulia/issues

In light of all the different issues, I've simply installed the command 
line version, added it to the Windows path (so I can type "julia" at any 
command prompt opened at any folder), and use the Atom editor to edit code. 
I have done this successfully on Windows 7 Julia 0.4 64-bit and I like this 
workflow because there have been zero issues.

I would recommend giving this a try to see if it fits your needs. An added 
benefit is that you can use any text editor you want to edit code.


On Sunday, November 1, 2015 at 11:31:31 PM UTC-5, John Wasa wrote:
>
> Donald,
>
> I installed Julia 4.0 and Sublime-IJulia on Windows 7 (64bit).
> Sublime shows that Ijulia kernel is working, but none of the shortcut keys 
> (shift+Enter, Ctrl+Enter, etc.) are working.
> Could you confirm that you also installed 64-bit Julia 4.0?  (My Version 
> is: 0.4.0 (2015-10-08 06:20 UTC) Official http://julialang.org/ release 
> x86_64-w64-mingw32)
>
> Thanks.
>


[julia-users] Re: [help]: configuring Julia and Sublime text 2 on windows 8.1

2015-10-03 Thread Donald Lacombe
Stefano,

I have had success using Sublime Text with Julia on Windows 7 using these 
instructions:

https://github.com/quinnj/Sublime-IJulia

If you follow these directions carefully, you should be able to get it up 
and running.

There are also two other alternatives which may appeal to you:

http://junolab.org/

https://github.com/JuliaLang/IJulia.jl

I'm not an expert but I have used all three options and they work just fine 
with the same basic functionality so the choice will boil down to personal 
preference.

On Saturday, October 3, 2015 at 11:17:57 AM UTC-4, stefano fedele wrote:
>
> Hello, I am completely new of Julia and mostly I am experienced with 
> interpreted programming languages. I would like to run Julia from Sublime 
> text 2 on Windows 8.1 but honestly I don't know where I should have to 
> start. Do I have to set the environment path for the Julia programming 
> language? Do I have to set something on the Sublime text 2 editor? Do I 
> have to attend other kind of procedures?
>
> Thank you in advance for your help :)
>


[julia-users] Re: Pretty Printing of Results

2014-09-24 Thread Donald Lacombe
Jason,

Thank you for this information. I will try to get this working soon. If you 
would like the code for the Bayesian spatial autoregressive model, please 
send me an email: donald.lacombe-at-mail.wvu.edu and I can get these and 
the support functions out to you.

The code currently uses simulated data but we can discuss the finalizing of 
the code if you are interested. I plan on coding other models as well.

Thanks,
Don



On Tuesday, September 23, 2014 8:45:12 PM UTC-4, Jason Knight wrote:

 Donald,

 I would enjoy seeing your code once you release it. I'm doing research 
 with Bayesian bioinformatics models using MCMC. 

 You can place all your MCMC samples in a DataFrame and then call describe 
 on it to get a very similar output to what you described.

 You can also write your own inspired from the definition 
 https://github.com/JuliaStats/DataFrames.jl/blob/cb9cd7f4c72eb49efd94d722731919a4086f683c/src/dataframe/dataframe.jl#L683
  
 of the describe function. Or instead of using rpad to get the spacing 
 right, you can use @printf or @sprintf instead.

 On Tuesday, September 23, 2014 6:29:25 PM UTC-5, Donald Lacombe wrote:

 Dear Julia Users,

 I have been coding various Bayesian spatial econometric models and have a 
 question that is probably very basic but the solution eludes me.

 I'd like to be able to pretty print the output from the models akin to 
 a standard regression package like so:

 Direct Effects
 VariablePosterior MeanLower 95%Upper 95%
 Income  3.0965  2.8765   3.4567
 Education   .3456 .2987 .3985

 I'd like to be able to have the user define the variables names and the 
 remaining quantities would be computed from the MCMC samplers.

 Can anyone provide a workable solution to this issue? Is there a general 
 purpose function where a user could input the variables names and matrix of 
 data to be printed?

 Any help that could be provided would be greatly appreciated. Also, if 
 anyone is interested in these models, please let me know. I'd like to 
 eventually release these for all to use after
 this issue is surmounted.

 Regards,
 Don 



Re: [julia-users] Re: Pretty Printing of Results

2014-09-24 Thread Donald Lacombe
Dear Milan,

Thank you for this information. I think that it will prove to be very 
useful for my needs.

Regards,
Don

On Wednesday, September 24, 2014 3:29:00 AM UTC-4, Milan Bouchet-Valat 
wrote:

  Le mardi 23 septembre 2014 à 17:45 -0700, Jason Knight a écrit : 

 Donald, 

  

  I would enjoy seeing your code once you release it. I'm doing research 
 with Bayesian bioinformatics models using MCMC.  

  

  You can place all your MCMC samples in a DataFrame and then call 
 describe on it to get a very similar output to what you described. 

  

  You can also write your own inspired from the definition 
 https://github.com/JuliaStats/DataFrames.jl/blob/cb9cd7f4c72eb49efd94d722731919a4086f683c/src/dataframe/dataframe.jl#L683
  
 of the describe function. Or instead of using rpad to get the spacing 
 right, you can use @printf or @sprintf instead. 

 You can also take inspiration from coeftable from StatsBase. One could 
 even imagine extending it to show either p-values (as is done currently) or 
 confidence intervals (as you apparently need), since different models 
 and/or different users are going to want to get either or both.

 Regards 



Re: [julia-users] Re: Pretty Printing of Results

2014-09-24 Thread Donald Lacombe
Dear Milan,

Do you have a link to the coeftable code? I looked here but could not 
find it:

https://github.com/JuliaStats/StatsBase.jl

Thanks,
Don

On Wednesday, September 24, 2014 3:29:00 AM UTC-4, Milan Bouchet-Valat 
wrote:

  Le mardi 23 septembre 2014 à 17:45 -0700, Jason Knight a écrit : 

 Donald, 

  

  I would enjoy seeing your code once you release it. I'm doing research 
 with Bayesian bioinformatics models using MCMC.  

  

  You can place all your MCMC samples in a DataFrame and then call 
 describe on it to get a very similar output to what you described. 

  

  You can also write your own inspired from the definition 
 https://github.com/JuliaStats/DataFrames.jl/blob/cb9cd7f4c72eb49efd94d722731919a4086f683c/src/dataframe/dataframe.jl#L683
  
 of the describe function. Or instead of using rpad to get the spacing 
 right, you can use @printf or @sprintf instead. 

 You can also take inspiration from coeftable from StatsBase. One could 
 even imagine extending it to show either p-values (as is done currently) or 
 confidence intervals (as you apparently need), since different models 
 and/or different users are going to want to get either or both.

 Regards 



Re: [julia-users] Re: Pretty Printing of Results

2014-09-24 Thread Donald Lacombe
Jason,

I copied lines 26-82 from the Git Hub page you kindly referenced and did 
the following:

In  [4]: mat = [1 2 3;4 5 6;7 8 9]

Out [4]: 3x3 Array{Int64,2}:
 1  2  3
 4  5  6
 7  8  9

In  [5]: colnms =[a;b;c]

Out [5]: 3-element Array{ASCIIString,1}:
 a
 b
 c

In  [6]: rownms = [a;b;c]

Out [6]: 3-element Array{ASCIIString,1}:
 a
 b
 c

In  [7]: CoefTable(mat,colnms,rownms)

Out [7]: CoefTable(3x3 Array{Int64,2}:
 1  2  3
 4  5  6
 7  8  9,ASCIIString[a,b,c],ASCIIString[a,b,c],0)

It seems to build the coefficient table through the CoefTable function but 
then I run into the following:

In  [11]: show(CoefTable)

Out [11]: ERROR: `show` has no method matching show(::Type{CoefTable})
you may have intended to import Base.show
while loading In[11], in expression starting on line 1


In  [12]: import Base.show

Warning: import of Base.show into Main conflicts with an existing 
identifier; ignored.

At this point I'm a bit stumped and was wondering if you or anyone else 
could advise me as to what to do. After I get this issue resolved, I think 
I should be able to alter the code to fit my needs.

Again, thank you for all of the help. I greatly appreciate the time and 
patience of everyone involved.

Regards,
Don


On Wednesday, September 24, 2014 2:53:45 PM UTC-4, Jason wrote:

 Do you have a link to the coeftable code? I looked here but could not 
 find it:


 Try: 
 https://github.com/JuliaStats/StatsBase.jl/blob/a2ba466f0b1ba82ed45f35918a9c5804ccc11d6b/src/statmodels.jl#L50
  



Re: [julia-users] Re: Pretty Printing of Results

2014-09-24 Thread Donald Lacombe
Jason,

I think you may have saved my sanity! Here is the output:

In  [22]: temp = CoefTable(mat,colnms,rownms)

Out [22]:  a b c
a1 2 3
b4 5 6
c7 8 9


In  [23]: show(temp)
 a b c
a1 2 3
b4 5 6
c7 8 9

These actually line up perfectly in my Sublime Editor so I think I'm fine. 
I also did the following as a more realistic example:

In  [30]: show(temp)
   Posterior Mean Lower 95% Upper 95%
Income  -0.451334   1.60543  -1.05135
Education 1.80312  0.382321  -2.06968
Price-0.58198  -1.53586  0.569491

Again, it lines up perfectly in my Sublime window. I think all I will have 
to do is to omit the p-value code and it should work fine with my spatial 
models.

Thank you again for answering my basic questions and this should make my 
code much nicer.

Regards,
Don

On Wednesday, September 24, 2014 8:22:23 PM UTC-4, Jason wrote:

 In  [11]: show(CoefTable)

 Out [11]: ERROR: `show` has no method matching show(::Type{CoefTable})
 you may have intended to import Base.show
 while loading In[11], in expression starting on line 1


 In  [12]: import Base.show

 Warning: import of Base.show into Main conflicts with an existing 
 identifier; ignored.


  Don, 

 You need to import Base.show *before* defining the show(::Type{CoefTable}) 
 that I'm assuming you copied from the link I sent earlier. Try doing it 
 that way and let us know how it goes.

 If you'd like to read more check out 
 http://docs.julialang.org/en/latest/manual/modules/



Re: [julia-users] Re: Pretty Printing of Results

2014-09-24 Thread Donald Lacombe
Dear Leah,

I actually used the following code that Jason told me about:

 
https://github.com/JuliaStats/StatsBase.jl/blob/a2ba466f0b1ba82ed45f35918a9c5804ccc11d6b/src/statmodels.jl#L50
 
https://www.google.com/url?q=https%3A%2F%2Fgithub.com%2FJuliaStats%2FStatsBase.jl%2Fblob%2Fa2ba466f0b1ba82ed45f35918a9c5804ccc11d6b%2Fsrc%2Fstatmodels.jl%23L50sa=Dsntz=1usg=AFQjCNEqvb02-hAiKV15uhONiZk2lTLaGQ
 

I post the code if you would like. I do not think I'm at the level of a 
JuliaExpert but I'd be happy to share what I did in my example to get 
others going.

Thoughts?

Regards,
Don

On Wednesday, September 24, 2014 8:52:28 PM UTC-4, Leah Hanson wrote:

 Hey Don,

 Have you considered putting your CoefTable function in a Julia package 
 somewhere?

 I'm not sure whether there's an existing one that would be appropriate, 
 but this looks like a useful bit of code that could go in a library.

 -- Leah

 On Wed, Sep 24, 2014 at 7:49 PM, Donald Lacombe drla...@gmail.com 
 javascript: wrote:

 Jason,

 I think you may have saved my sanity! Here is the output:

 In  [22]: temp = CoefTable(mat,colnms,rownms)

 Out [22]:  a b c
 a1 2 3
 b4 5 6
 c7 8 9


 In  [23]: show(temp)
  a b c
 a1 2 3
 b4 5 6
 c7 8 9

 These actually line up perfectly in my Sublime Editor so I think I'm 
 fine. I also did the following as a more realistic example:

 In  [30]: show(temp)
Posterior Mean Lower 95% Upper 95%
 Income  -0.451334   1.60543  -1.05135
 Education 1.80312  0.382321  -2.06968
 Price-0.58198  -1.53586  0.569491

 Again, it lines up perfectly in my Sublime window. I think all I will 
 have to do is to omit the p-value code and it should work fine with my 
 spatial models.

 Thank you again for answering my basic questions and this should make my 
 code much nicer.

 Regards,
 Don

 On Wednesday, September 24, 2014 8:22:23 PM UTC-4, Jason wrote:

 In  [11]: show(CoefTable)

 Out [11]: ERROR: `show` has no method matching show(::Type{CoefTable})
 you may have intended to import Base.show
 while loading In[11], in expression starting on line 1


 In  [12]: import Base.show

 Warning: import of Base.show into Main conflicts with an existing 
 identifier; ignored.


  Don, 

 You need to import Base.show *before* defining the 
 show(::Type{CoefTable}) that I'm assuming you copied from the link I sent 
 earlier. Try doing it that way and let us know how it goes.

 If you'd like to read more check out http://docs.julialang.org/
 en/latest/manual/modules/




[julia-users] Pretty Printing of Results

2014-09-23 Thread Donald Lacombe
Dear Julia Users,

I have been coding various Bayesian spatial econometric models and have a 
question that is probably very basic but the solution eludes me.

I'd like to be able to pretty print the output from the models akin to a 
standard regression package like so:

Direct Effects
VariablePosterior MeanLower 95%Upper 95%
Income  3.0965  2.8765   3.4567
Education   .3456 .2987 .3985

I'd like to be able to have the user define the variables names and the 
remaining quantities would be computed from the MCMC samplers.

Can anyone provide a workable solution to this issue? Is there a general 
purpose function where a user could input the variables names and matrix of 
data to be printed?

Any help that could be provided would be greatly appreciated. Also, if 
anyone is interested in these models, please let me know. I'd like to 
eventually release these for all to use after
this issue is surmounted.

Regards,
Don 


[julia-users] Re: Julia equivalent of MATLAB colamd

2014-08-01 Thread Donald Lacombe
Dear Doug,

My apologies for misstating something. The variance-covariance matrix for 
the errors in these spatial models is positive definite, not the (In - 
rho*W) term.

The W matrix is not symmetric in general. My goal is to create a lookup 
table of values of log(det(In-rho*W)) for values of rho from lmin to lmax 
and then simply look these values up during the Metropolis step. This is 
why I asked if the code could be used in the case of a non-symmetric W 
matrix.

Thank you again for taking the time to look at this.

Don



On Friday, August 1, 2014 3:20:11 PM UTC-4, Douglas Bates wrote:

 On Thursday, July 31, 2014 8:33:49 PM UTC-5, Donald Lacombe wrote:

 Dear Doug,

 Thank you for taking the time to look at this. In actuality, the W matrix 
 is a spatial weight matrix which defines who is a neighbor to whom. I think 
 the symmetric verbiage in the code is not technically correct. The 
 standardized means that  the rows of the W matrix sum to 1. This 
 ensures that the (In -rho*W) matrix is positive definite as you state.


 I'm confused about the terminology then.  To me it doesn't make sense to 
 say that (I-rho*W) is positive definite unless W is symmetric (or 
 Hermitian, for complex arithmetic) as stated in 
 http://en.wikipedia.org/wiki/Positive-definite_matrix.
  


 The following matrix would not technically be a proper spatial weight 
 matrix because the 4th row has 1 neighbor and the 5th row has three 
 neighbors when in fact they should all have 2 neighbors (actually there can 
 be a different number of neighbors but I'm choosing 2 just for an 
 illustration).

 julia full(W)
 5x5 Array{Float64,2}:
  0.0  0.5  0.5  0.0  0.0
  0.5  0.0  0.0  0.0  0.5
  0.5  0.0  0.0  0.0  0.5
  0.0  0.0  0.0  0.0  0.5
  0.0  0.5  0.5  0.5  0.0

 The original matrix I posted was a proper nearest neighbor spatial weight 
 matrix:

  full(W)

 ans =

  00.50000.5000 0 0
 0.5000 0 0 00.5000
 0.5000 0 0 00.5000
  00.5000 0 00.5000
  0 00.50000.5000 0 

 Each row represent a geographic entity and the columns represent their 
 neighbors. The diagonal elements are zero because no geographic entity is a 
 neighbor to itself. The 0.5's represent the row-normalization and in this 
 case the value is .5 because there are 2 neighbors.

 I guess the bottom line is that the code you provided works but I would 
 need it to work for a non-symmetric W matrix.

 Can the above code be altered to deal with the non-symmetric W case?

 Thank you again for taking the time to look at this. I truly appreciate 
 the help and education.

 Regards,
 Don


 On Thursday, July 31, 2014 4:49:21 PM UTC-4, Douglas Bates wrote:

 Rats, I had a reply with code and examples then Google groups croaked 
 and when it reloaded my draft reply was gone.

 On Thursday, July 31, 2014 2:13:36 PM UTC-5, Donald Lacombe wrote:

 Dear Doug,

 I've created a small example to show what the MATLAB code is actually 
 doing. Here is the output with some comments:

 % Create random coordinates
  xc = randn(5,1);
  yc = randn(5,1);

 % Create 2 nearest neighbor weight matrix
  W = make_neighborsw(xc,yc,2);
  W

 W =

(2,1)   0.5000
(3,1)   0.5000
(1,2)   0.5000
(4,2)   0.5000
(1,3)   0.5000
(5,3)   0.5000
(5,4)   0.5000
(2,5)   0.5000
(3,5)   0.5000
(4,5)   0.5000

 % Full version
  full(W)

 ans =

  00.50000.5000 0 0
 0.5000 0 0 00.5000
 0.5000 0 0 00.5000
  00.5000 0 00.5000
  0 00.50000.5000 0


 Notice that this matrix is not symmetric but the description of the 
 arguments for indetfull below says that W should be symmetric.  I am going 
 to assume that W should be symmetric and the allowable values of rho are 
 such that I - rho*W is positive-definite. I will use the matrix  

 julia W = sparse([2,3,1,5,1,5,5,2,3,4],[1,1,2,2,3,3,4,5,5,5],0.5)
 5x5 sparse matrix with 10 Float64 entries:
 [2, 1]  =  0.5
 [3, 1]  =  0.5
 [1, 2]  =  0.5
 [5, 2]  =  0.5
 [1, 3]  =  0.5
 [5, 3]  =  0.5
 [5, 4]  =  0.5
 [2, 5]  =  0.5
 [3, 5]  =  0.5
 [4, 5]  =  0.5

 julia issym(W)
 true

 julia full(W)
 5x5 Array{Float64,2}:
  0.0  0.5  0.5  0.0  0.0
  0.5  0.0  0.0  0.0  0.5
  0.5  0.0  0.0  0.0  0.5
  0.0  0.0  0.0  0.0  0.5
  0.0  0.5  0.5  0.5  0.0


 for illustration.


 % Find colamd
  z = speye(n) - 0.1*sparse(W)

 z =

(1,1)   1.
(2,1)  -0.0500
(3,1)  -0.0500
(1,2)  -0.0500
(2,2)   1.
(4,2)  -0.0500
(1,3)  -0.0500
(3,3)   1.
(5,3)  -0.0500
(4,4)   1.
(5,4)  -0.0500
(2,5)  -0.0500
(3,5)  -0.0500
(4,5)  -0.0500
(5,5)   1.

 % Full z
  full

[julia-users] Re: Julia equivalent of MATLAB colamd

2014-07-31 Thread Donald Lacombe
Dear Doug,

I've created a small example to show what the MATLAB code is actually 
doing. Here is the output with some comments:

% Create random coordinates
 xc = randn(5,1);
 yc = randn(5,1);

% Create 2 nearest neighbor weight matrix
 W = make_neighborsw(xc,yc,2);
 W

W =

   (2,1)   0.5000
   (3,1)   0.5000
   (1,2)   0.5000
   (4,2)   0.5000
   (1,3)   0.5000
   (5,3)   0.5000
   (5,4)   0.5000
   (2,5)   0.5000
   (3,5)   0.5000
   (4,5)   0.5000

% Full version
 full(W)

ans =

 00.50000.5000 0 0
0.5000 0 0 00.5000
0.5000 0 0 00.5000
 00.5000 0 00.5000
 0 00.50000.5000 0

% Find colamd
 z = speye(n) - 0.1*sparse(W)

z =

   (1,1)   1.
   (2,1)  -0.0500
   (3,1)  -0.0500
   (1,2)  -0.0500
   (2,2)   1.
   (4,2)  -0.0500
   (1,3)  -0.0500
   (3,3)   1.
   (5,3)  -0.0500
   (4,4)   1.
   (5,4)  -0.0500
   (2,5)  -0.0500
   (3,5)  -0.0500
   (4,5)  -0.0500
   (5,5)   1.

% Full z
 full(z)

ans =

1.   -0.0500   -0.0500 0 0
   -0.05001. 0 0   -0.0500
   -0.0500 01. 0   -0.0500
 0   -0.0500 01.   -0.0500
 0 0   -0.0500   -0.05001.

 p = colamd(z)

p =

 4 3 1 2 5

What I'm looking for is the p vector above. The full MATLAB code for the 
routine is as follows:

function out=lndetfull(W,lmin,lmax)
% PURPOSE: computes Pace and Barry's grid for log det(I-rho*W) using sparse 
matrices
% ---
% USAGE: out = lndetfull(W,lmin,lmax)
% where:
% W = symmetric spatial weight matrix (standardized)
% lmin  = lower bound on rho
% lmax  = upper bound on rho
% ---
% RETURNS: out = a structure variable
%  out.lndet = a vector of log determinants for 0  rho  1
%  out.rho   = a vector of rho values associated with lndet values
% ---
% NOTES: should use 1/lambda(max) to 1/lambda(min) for all possible rho 
values
% ---
% References: % R. Kelley Pace and  Ronald Barry. 1997. ``Quick
% Computation of Spatial Autoregressive Estimators'', Geographical Analysis
% ---
 
% written by:
% James P. LeSage, Dept of Economics
% University of Toledo
% 2801 W. Bancroft St,
% Toledo, OH 43606


rvec = lmin:.01:lmax;
spparms('tight'); 
[n junk] = size(W);
z = speye(n) - 0.1*sparse(W);
p = colamd(z);
niter = length(rvec);
dettmp = zeros(niter,2);
for i=1:niter;
rho = rvec(i);
z = speye(n) - rho*sparse(W);
[l,u] = lu(z(:,p));
dettmp(i,1) = rho;
dettmp(i,2) = sum(log(abs(diag(u;
end;

out.lndet = dettmp(:,2);
out.rho = dettmp(:,1);

The code will be very easy to convert once I figure out this one command,

Thank you in advance for any help that you are able or willing to provide. 
I'm not an expert as sparse matrices but these routines come in very handy 
for these types of calculations.

Thanks,
Don


On Thursday, July 31, 2014 10:42:55 AM UTC-4, Douglas Bates wrote:

 On Wednesday, July 30, 2014 8:10:37 PM UTC-5, Donald Lacombe wrote:

 Dear Julia Users:

 I am attempting to speed up some code regarding some Bayesian spatial 
 econometric models and would like to translate some MATLAB code for 
 computing the log determinant term in the log-likelihood function. 
 Specifically, I'm interested in the Pace and Barry approximation to the 
 log-determinant term as in R. Kelley Pace and Ronald Barry. (1997) Quick 
 Computation 
 of Spatial Autoregressive Estimators''. The term in question is det(In - 
 rho*W) which has to be calculated through each iteration of the 
 Metropolis-Hastings algorithm I'm using.

 The MATLAB code for this has the command colamd which MATLAB's help 
 defines as P = colamd(S) returns the column approximate minimum degree 
 permutation vector for the sparse matrix S.


 colamd is indirectly available in the Cholmod module but I don't think 
 there is a direct way of calling it.  It is a method for determining a 
 fill-reducing permutation for a sparse, symmetric matrix of the form X'X 
 directly from X.  There are several other approaches implemented in amd, 
 Metis, Scotch, MUMPS, etc. for dealing with the symmetric matrix itself. 
  All of these just use the pattern of nonzeros in the matrix as they are 
 part of the symbolic stage of the sparse Cholesky or sparse QR 
 factorization. 

 Can you be more specific about how it is to be used?   In particular, is 
 it important to work

[julia-users] Julia equivalent of MATLAB colamd

2014-07-30 Thread Donald Lacombe
Dear Julia Users:

I am attempting to speed up some code regarding some Bayesian spatial 
econometric models and would like to translate some MATLAB code for 
computing the log determinant term in the log-likelihood function. 
Specifically, I'm interested in the Pace and Barry approximation to the 
log-determinant term as in R. Kelley Pace and Ronald Barry. (1997) Quick 
Computation 
of Spatial Autoregressive Estimators''. The term in question is det(In - 
rho*W) which has to be calculated through each iteration of the 
Metropolis-Hastings algorithm I'm using.

The MATLAB code for this has the command colamd which MATLAB's help 
defines as P = colamd(S) returns the column approximate minimum degree 
permutation vector for the sparse matrix S.

Is there an equivalent Julia command? I've tried looking at symperm but I 
do not think I am using it correctly and more importantly I have no real 
clue.

Any help would be appreciated and if I successfully code this I would be 
happy to share with anyone who is interested.

Thanks,
Don 


[julia-users] Question Regarding Issue Tracker

2014-07-26 Thread Donald Lacombe
Dear Julia Users,

A while ago, Simon Byrne files the following issue on GitHub:

https://github.com/JuliaLang/julia/issues/7543

I see that the issue is currently closed but I'm unsure if the suggestion 
has been incorporated or is just there for posterity.

Is there a way to find out what happened? If this is a newbie GitHub 
issue my apologies in advance. If there is a tutorial that addresses this 
it might be worth posting a link.

Thank you in advance for any help you can provide.

Regards,
Don


Re: [julia-users] Question Regarding Issue Tracker

2014-07-26 Thread Donald Lacombe
Dear John,

Thank you for that bit of information! 

One further query: I assume that these fixes appear in the daily builds and 
will be incorporated into the next Current Release?

Thank you agian for the quick response and the help.

Regards,
Don

On Saturday, July 26, 2014 12:10:20 PM UTC-4, John Myles White wrote:

 If you look carefully at the issue where the “merged” icon occurs, you’ll 
 see that a patch was merged that solved the issue. This is a new way that 
 GitHub has recently started visualizing events inside of an issue comment 
 thread. 

  — John 

 On Jul 26, 2014, at 9:07 AM, Donald Lacombe drla...@gmail.com 
 javascript: wrote: 

  Dear Julia Users, 
  
  A while ago, Simon Byrne files the following issue on GitHub: 
  
  https://github.com/JuliaLang/julia/issues/7543 
  
  I see that the issue is currently closed but I'm unsure if the 
 suggestion has been incorporated or is just there for posterity. 
  
  Is there a way to find out what happened? If this is a newbie GitHub 
 issue my apologies in advance. If there is a tutorial that addresses this 
 it might be worth posting a link. 
  
  Thank you in advance for any help you can provide. 
  
  Regards, 
  Don 



[julia-users] Re: Using Sparse Matrix Algorithms: Spatial Econometrics

2014-07-11 Thread Donald Lacombe
Simon,

Thanks again for the input. I searched the julia stats forum but there 
did not seem to be any posts regarding spatial econometric models.

I might post the revised code there to see if there is any mutual interest. 
I have a few more items to code and the final piece will be a pretty 
printing function somewhat akin to the output of R, etc.

Thanks,
Don

On Wednesday, July 9, 2014 5:23:19 AM UTC-4, Simon Byrne wrote:

 Hi Don,

 On Tuesday, 8 July 2014 20:51:42 UTC+1, Donald Lacombe wrote:

 2) When I use A = sparse(In - rho*W) I get the following warning from 
 Julia and a break:

 julia sar_gibbs

 LoadError(C:/Users/Don/Desktop/Julia Materials/SAR 
 Gibbs/sar_gibbs.jl,13,ErrorException(The inverse of a sparse matrix can 
 often be dense and can cause the computer to run out of memory. If you are 
 sure you have enough memory, please convert your matrix to a dense matrix.))


 Apparently, it does not like it when I invert the A matrix.


 For reasons of performance and numerical accuracy, you should generally 
 try to avoid explicitly constructing matrix inverses where possible: 
 instead of writing inv(A) * y, you should use the left division operator A 
 \ y.

 3) Can I wrap the code with a blank function call like so: function 
 sar_gibbs()? 


 Yes, that's what I meant.
  

 4) Are you indicating that I can include the knn_weights_matrix code in the 
 same file?


 Yes: obviously this is up to you, but for short functions like this I find 
 it much easier to keep track of everything in the same file.
  

 Once again, thank you for the response. My goal is to produce a suite of 
 Bayesian spatial econometric models for public consumption and you help is 
 greatly appreciated.


 There are quite a few people interested in similar models: you might also 
 be interested in subscribing to the julia-stats mailing list.

 Good luck,

 Simon

  

 On Tuesday, July 8, 2014 5:06:14 AM UTC-4, Simon Byrne wrote:

 Hi Don,

 I think the reason they're not sparse is that in the line
 A = (In - rho*W);
 the matrix In is not sparse: if you replace the line
 In = eye(n);
 by 
 In = speye(n);
 the result should then be sparse. However at the moment there doesn't 
 seem to be a sparse det method (I just filed issue #7543 
 https://github.com/JuliaLang/julia/issues/7543): you can get around 
 this by using log(det(lufact(A))) instead of log(det(A)).

 You might also want to have a quick look at the performance tips 
 http://docs.julialang.org/en/latest/manual/performance-tips/, in 
 particular:

 * You should use randn() instead of randn(1): randn() samples a single 
 Float64 scalar, randn(1) samples a vector of length 1. On my machine (using 
 0.3 pre-release) your code doesn't work, as this promotes some vectors to 
 matrices, which doesn't work with dot. You also then won't need to call 
 variables like rho_proposed[1].
 * Try wrapping everything in a function: this makes it easier for the 
 JIT compiler to work its magic, as well as making it easier to use 
 profiling tools.
 * If you're using dense matrices, you can use logdet(A) instead of 
 log(det(A)): this can avoid some underflow and overflow problems 
 (unfortunately this doesn't work for sparse lu factorizations yet, see the 
 issue mentioned above).
 * Purely from a style point of view, there's no need to keep functions 
 in separate files, or end lines with semi-colons.

 Simon



 On Tuesday, 8 July 2014 01:48:56 UTC+1, Donald Lacombe wrote:

 Dear Julia Users,

 I am currently developing/converting several Bayesian spatial 
 econometric models from MATLAB into Julia. I have successfully coded the 
 spatial autoregressive model (SAR) with diffuse priors in Julia but have a 
 question regarding the use of sparse matrix algorithms. The models use a 
 spatial weight matrix which is usually sparse in nature and this matrix 
 appears in many of the expressions, especially in the random walk 
 Metropolis algorithm used to obtain the marginal distribution for the 
 spatial autocorrelation parameter rho. Here is a code snippet and the 
 complete code is attached:

 # Draw for rho

 A = (In - rho*W);

 denominator = log(det(A))-.5*sigma^-2*dot(A*y-x*beta,A*y-x*beta);

 accept = 0;

 rho_proposed = rho + zta*randn(1);

 while accept == 0

 if ((rho_proposed[1]  -1)  (rho_proposed[1]  1));

 accept = 1;

 else

 rho_proposed = rho + zta*randn(1);

 end

 end

  B = (In - rho_proposed[1]*W);

  numerator = log(det(B))-.5*sigma^-2*dot(B*y-x*beta,B*y-x*beta);

  u = rand(1);

  if ((numerator[1] - denominator[1])  exp(1))

 pp = 1;

  else

 ratio = exp(numerator[1] - denominator[1]);

 pp = min(1,ratio);

 end

 if (u[1]  pp[1])

 rho = rho_proposed[1];

 acc = acc + 1;

 end

 ar = acc/gibbs;


 if ar  0.4

 zta = zta/1.1;

 end

 if ar  0.6

 zta = zta*1.1;

 end

[julia-users] Re: Using Sparse Matrix Algorithms: Spatial Econometrics

2014-07-08 Thread Donald Lacombe
Dear Simon,

Thank you for the quick and informative response! I have a few 
questions/observations:

1) Thank you for the information regarding the randn() function. I found a 
workaround in the original code but now it looks much better and acts as 
expected.

2) When I use A = sparse(In - rho*W) I get the following warning from Julia 
and a break:

julia sar_gibbs

LoadError(C:/Users/Don/Desktop/Julia Materials/SAR 
Gibbs/sar_gibbs.jl,13,ErrorException(The inverse of a sparse matrix can often 
be dense and can cause the computer to run out of memory. If you are sure you 
have enough memory, please convert your matrix to a dense matrix.))


Apparently, it does not like it when I invert the A matrix.


3) Can I wrap the code with a blank function call like so: function 
sar_gibbs()? 


4) Are you indicating that I can include the knn_weights_matrix code in the 
same file?


Once again, thank you for the response. My goal is to produce a suite of 
Bayesian spatial econometric models for public consumption and you help is 
greatly appreciated.


Regards,

Don




On Tuesday, July 8, 2014 5:06:14 AM UTC-4, Simon Byrne wrote:

 Hi Don,

 I think the reason they're not sparse is that in the line
 A = (In - rho*W);
 the matrix In is not sparse: if you replace the line
 In = eye(n);
 by 
 In = speye(n);
 the result should then be sparse. However at the moment there doesn't seem 
 to be a sparse det method (I just filed issue #7543 
 https://github.com/JuliaLang/julia/issues/7543): you can get around 
 this by using log(det(lufact(A))) instead of log(det(A)).

 You might also want to have a quick look at the performance tips 
 http://docs.julialang.org/en/latest/manual/performance-tips/, in 
 particular:

 * You should use randn() instead of randn(1): randn() samples a single 
 Float64 scalar, randn(1) samples a vector of length 1. On my machine (using 
 0.3 pre-release) your code doesn't work, as this promotes some vectors to 
 matrices, which doesn't work with dot. You also then won't need to call 
 variables like rho_proposed[1].
 * Try wrapping everything in a function: this makes it easier for the JIT 
 compiler to work its magic, as well as making it easier to use profiling 
 tools.
 * If you're using dense matrices, you can use logdet(A) instead of 
 log(det(A)): this can avoid some underflow and overflow problems 
 (unfortunately this doesn't work for sparse lu factorizations yet, see the 
 issue mentioned above).
 * Purely from a style point of view, there's no need to keep functions in 
 separate files, or end lines with semi-colons.

 Simon



 On Tuesday, 8 July 2014 01:48:56 UTC+1, Donald Lacombe wrote:

 Dear Julia Users,

 I am currently developing/converting several Bayesian spatial econometric 
 models from MATLAB into Julia. I have successfully coded the spatial 
 autoregressive model (SAR) with diffuse priors in Julia but have a question 
 regarding the use of sparse matrix algorithms. The models use a spatial 
 weight matrix which is usually sparse in nature and this matrix appears in 
 many of the expressions, especially in the random walk Metropolis algorithm 
 used to obtain the marginal distribution for the spatial autocorrelation 
 parameter rho. Here is a code snippet and the complete code is attached:

 # Draw for rho

 A = (In - rho*W);

 denominator = log(det(A))-.5*sigma^-2*dot(A*y-x*beta,A*y-x*beta);

 accept = 0;

 rho_proposed = rho + zta*randn(1);

 while accept == 0

 if ((rho_proposed[1]  -1)  (rho_proposed[1]  1));

 accept = 1;

 else

 rho_proposed = rho + zta*randn(1);

 end

 end

  B = (In - rho_proposed[1]*W);

  numerator = log(det(B))-.5*sigma^-2*dot(B*y-x*beta,B*y-x*beta);

  u = rand(1);

  if ((numerator[1] - denominator[1])  exp(1))

 pp = 1;

  else

 ratio = exp(numerator[1] - denominator[1]);

 pp = min(1,ratio);

 end

 if (u[1]  pp[1])

 rho = rho_proposed[1];

 acc = acc + 1;

 end

 ar = acc/gibbs;


 if ar  0.4

 zta = zta/1.1;

 end

 if ar  0.6

 zta = zta*1.1;

 end


 A = (In - rho*W);


 Basically, I'm wondering if there is any way to make the A and B 
 matrices sparse and possibly make it run faster, especially in the 
 log(det(A)) terms. Currently, 110 draws (with 10 burn) takes approximately 8 
 seconds on my 64 bit, core i7 laptop. The computational speed decreases with 
 the sample size n because the weight matrices are treated as full.


 Any help would be greatly appreciated and if anyone is interested in running 
 the code, the Distributions and Distance packages must be included and 
 initiated first.


 Regards,

 Don





Re: [julia-users] Re: Using Sparse Matrix Algorithms: Spatial Econometrics

2014-07-08 Thread Donald Lacombe
Andreas,

Yes, you are correct and your solution works just fine. 

Thank you for taking the time to answer my question. I've said it before 
but the kindness and generosity of this community is remarkable!

Thanks,
Don

On Tuesday, July 8, 2014 5:34:45 AM UTC-4, Andreas Noack wrote:

 I think I-ρW would be even better in this case. I is a generic identity 
 matrix. A and B are probably positive definite. If so, I think 
 logdet(cholfact(A)) would be the best solution. 


 2014-07-08 11:06 GMT+02:00 Simon Byrne simon...@gmail.com javascript::

 Hi Don,

 I think the reason they're not sparse is that in the line
 A = (In - rho*W);
 the matrix In is not sparse: if you replace the line
 In = eye(n);
 by 
 In = speye(n);
 the result should then be sparse. However at the moment there doesn't 
 seem to be a sparse det method (I just filed issue #7543 
 https://github.com/JuliaLang/julia/issues/7543): you can get around 
 this by using log(det(lufact(A))) instead of log(det(A)).

 You might also want to have a quick look at the performance tips 
 http://docs.julialang.org/en/latest/manual/performance-tips/, in 
 particular:

 * You should use randn() instead of randn(1): randn() samples a single 
 Float64 scalar, randn(1) samples a vector of length 1. On my machine (using 
 0.3 pre-release) your code doesn't work, as this promotes some vectors to 
 matrices, which doesn't work with dot. You also then won't need to call 
 variables like rho_proposed[1].
 * Try wrapping everything in a function: this makes it easier for the JIT 
 compiler to work its magic, as well as making it easier to use profiling 
 tools.
 * If you're using dense matrices, you can use logdet(A) instead of 
 log(det(A)): this can avoid some underflow and overflow problems 
 (unfortunately this doesn't work for sparse lu factorizations yet, see the 
 issue mentioned above).
 * Purely from a style point of view, there's no need to keep functions in 
 separate files, or end lines with semi-colons.

 Simon



 On Tuesday, 8 July 2014 01:48:56 UTC+1, Donald Lacombe wrote:

 Dear Julia Users,

 I am currently developing/converting several Bayesian spatial 
 econometric models from MATLAB into Julia. I have successfully coded the 
 spatial autoregressive model (SAR) with diffuse priors in Julia but have a 
 question regarding the use of sparse matrix algorithms. The models use a 
 spatial weight matrix which is usually sparse in nature and this matrix 
 appears in many of the expressions, especially in the random walk 
 Metropolis algorithm used to obtain the marginal distribution for the 
 spatial autocorrelation parameter rho. Here is a code snippet and the 
 complete code is attached:

 # Draw for rho

 A = (In - rho*W);

 denominator = log(det(A))-.5*sigma^-2*dot(A*y-x*beta,A*y-x*beta);

 accept = 0;

 rho_proposed = rho + zta*randn(1);

 while accept == 0

 if ((rho_proposed[1]  -1)  (rho_proposed[1]  1));

 accept = 1;

 else

 rho_proposed = rho + zta*randn(1);

 end

 end

  B = (In - rho_proposed[1]*W);

  numerator = log(det(B))-.5*sigma^-2*dot(B*y-x*beta,B*y-x*beta);

  u = rand(1);

  if ((numerator[1] - denominator[1])  exp(1))

 pp = 1;

  else

 ratio = exp(numerator[1] - denominator[1]);

 pp = min(1,ratio);

 end

 if (u[1]  pp[1])

 rho = rho_proposed[1];

 acc = acc + 1;

 end

 ar = acc/gibbs;


 if ar  0.4

 zta = zta/1.1;

 end

 if ar  0.6

 zta = zta*1.1;

 end


 A = (In - rho*W);


 Basically, I'm wondering if there is any way to make the A and B 
 matrices sparse and possibly make it run faster, especially in the 
 log(det(A)) terms. Currently, 110 draws (with 10 burn) takes approximately 
 8 seconds on my 64 bit, core i7 laptop. The computational speed decreases 
 with the sample size n because the weight matrices are treated as full.


 Any help would be greatly appreciated and if anyone is interested in 
 running the code, the Distributions and Distance packages must be included 
 and initiated first.


 Regards,

 Don






 -- 
 Med venlig hilsen

 Andreas Noack Jensen
  


[julia-users] Using Sparse Matrix Algorithms: Spatial Econometrics

2014-07-07 Thread Donald Lacombe
Dear Julia Users,

I am currently developing/converting several Bayesian spatial econometric 
models from MATLAB into Julia. I have successfully coded the spatial 
autoregressive model (SAR) with diffuse priors in Julia but have a question 
regarding the use of sparse matrix algorithms. The models use a spatial 
weight matrix which is usually sparse in nature and this matrix appears in 
many of the expressions, especially in the random walk Metropolis algorithm 
used to obtain the marginal distribution for the spatial autocorrelation 
parameter rho. Here is a code snippet and the complete code is attached:

# Draw for rho

A = (In - rho*W);

denominator = log(det(A))-.5*sigma^-2*dot(A*y-x*beta,A*y-x*beta);

accept = 0;

rho_proposed = rho + zta*randn(1);

while accept == 0

if ((rho_proposed[1]  -1)  (rho_proposed[1]  1));

accept = 1;

else

rho_proposed = rho + zta*randn(1);

end

end

 B = (In - rho_proposed[1]*W);

 numerator = log(det(B))-.5*sigma^-2*dot(B*y-x*beta,B*y-x*beta);

 u = rand(1);

 if ((numerator[1] - denominator[1])  exp(1))

pp = 1;

 else

ratio = exp(numerator[1] - denominator[1]);

pp = min(1,ratio);

end

if (u[1]  pp[1])

rho = rho_proposed[1];

acc = acc + 1;

end

ar = acc/gibbs;


if ar  0.4

zta = zta/1.1;

end

if ar  0.6

zta = zta*1.1;

end


A = (In - rho*W);


Basically, I'm wondering if there is any way to make the A and B matrices 
sparse and possibly make it run faster, especially in the log(det(A)) terms. 
Currently, 110 draws (with 10 burn) takes approximately 8 seconds on my 64 bit, 
core i7 laptop. The computational speed decreases with the sample size n 
because the weight matrices are treated as full.


Any help would be greatly appreciated and if anyone is interested in running 
the code, the Distributions and Distance packages must be included and 
initiated first.


Regards,

Don





sar_gibbs.jl
Description: Binary data


knn_weight_matrix.jl
Description: Binary data


[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-03 Thread Donald Lacombe
Dear Tony,

Thank you for pointing this out. Worked like a charm and the distance were 
calculated just fine.

Thanks,
Don

On Wednesday, July 2, 2014 9:46:06 PM UTC-4, Tony Kelman wrote:

 If you know all of the data in the csv file is numerical, you can provide 
 a result type to readcsv, say readcsv(ct_coord.csv, Float64). See 
 help(readdlm) for additional options for things like removing header rows, 
 etc. If the data is of mixed types, then you may find readtable in the 
 DataFrames package more capable.


 On Wednesday, July 2, 2014 5:05:21 PM UTC-7, Donald Lacombe wrote:

 Actually, another issue with the type information.

 Instead of using random coordinates, I used actual latitude and longitude 
 coordinates from a CSV file like so:

 data = readcsv(ct_coord.csv);

 temp = [xc';yc']

 2x8 Array{Any,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


 Now, the type here is {Any,2} and this seems to be giving the Distance 
 package problems:



 R = pairwise(Euclidean(),temp)

 MethodError(fptype,(Any,))


 Julia seems to be unhappy with this type but I'm not sure how to correct 
 this which would seem to require that the data being read in are of Float 
 type.


 Again, thank you for all of the help.


 Don



 On Wednesday, July 2, 2014 7:10:09 PM UTC-4, Donald Lacombe wrote:

 Patrick,

 Thank you for the reply! I will try your suggestion and hopefully get a 
 better grasp of these issues.

 Thank you all gain for all of the help. It cannot be stressed enough 
 that the friendly and truly helpful comments will help the Julia community 
 grow.

 Don

 On Wednesday, July 2, 2014 5:31:55 PM UTC-4, Patrick O'Leary wrote:

 The second parameter to the Array type is the number of dimensions of 
 the array. Row vectors are implemented with 2-dimensional arrays 
 (a.k.a., 
 matrices), since Vectors (a.k.a. Array{T, 1} for all types T) are 
 interpreted as columns.

 Since pairwise() generates all of the distances simultaneously, you'll 
 need to slice its results to use it in sortperm(), following your original 
 loop implementation. Each slice will need to be a Vector. If you slice 
 columns, you get Vector slices for free. If you slice rows, the shape will 
 be maintained, so you get a row back. You can flatten this to a Vector 
 with 
 the vec() function.

 Hope this helps!

 On Wednesday, July 2, 2014 4:23:52 PM UTC-5, Donald Lacombe wrote:

 After thinking about this issue, another couple of question came to 
 mind.

 When I create a row vector, say xc = randn(1,5), the results indicate 
 that the vector is e.g.

 xc = rand(1,5)

 1x5 Array{Float64,2}:

  0.989421  0.718148  0.499914  0.13024  0.939856


 The type is {Float64,2}.


 Questions:


 1) What does the 2 mean in {Float64,2}? If I create xc = rand(5) the 
 type is {Float64,1}.


 2) The sortperm function only seems to accept {Float64,1} types. The 
 Euclidean distance code create types of {Float64,2} but I need 
 {Float64,1} for the sortperm function. Do I need to convert these before 
 using this command?


 I'm probably missing something really basic so apologies for the query.


 Thanks,

 Don


 On Wednesday, July 2, 2014 4:20:55 PM UTC-4, Donald Lacombe wrote:

 Greetings, part 2...

 After doing some tests, I found that this code works just fine:

 julia using Distance

 Warning: could not import Base.foldl into NumericExtensions

 Warning: could not import Base.foldr into NumericExtensions

 Warning: could not import Base.sum! into NumericExtensions

 Warning: could not import Base.maximum! into NumericExtensions

 Warning: could not import Base.minimum! into NumericExtensions



 julia xc = rand(8)'

 1x8 Array{Float64,2}:

  0.30107  0.479169  0.230607  0.65126  0.386921  0.455316  0.921496  
 0.244873


 julia yc = rand(8)'

 1x8 Array{Float64,2}:

  0.866199  0.767794  0.76163  0.262925  …  0.742765  0.980952  0.424966


 julia temp = [xc ; yc]

 2x8 Array{Float64,2}:

  0.30107   0.479169  0.230607  0.65126   …  0.455316  0.921496  0.244873

  0.866199  0.767794  0.76163   0.262925 0.742765  0.980952  0.424966


 julia R = pairwise(Euclidean(),temp)

 8x8 Array{Float64,2}:

  0.0   0.203477   0.126094  0.697548  …  0.197554   0.630949  
 0.444797

  0.203477  0.00.248638  0.533393 0.0345751  0.491009  
 0.415241

  0.126094  0.248638   0.0   0.652423 0.225499   0.724865  
 0.336966

  0.697548  0.533393   0.652423  0.0  0.518306   0.767196  
 0.437502

  0.376201  0.283308   0.304834  0.355027 0.252287   0.719136  
 0.160613

  0.197554  0.0345751  0.225499  0.518306  …  0.00.523505  
 0.381159

  0.630949  0.491009   0.724865  0.767196 0.523505   0.0   
 0.87575 

  0.444797  0.415241   0.336966  0.437502 0.381159   0.87575   0.0
  


 As expected, the diagonal elements are 0 and everything seems to work 
 fine although I did not test other distance formulas

[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-03 Thread Donald Lacombe
Dear Patrick,

Yes, I agree...I actually have some MATLAB code that calculates pairwise 
distances using a Great Circle formula. The Econometrics Toolbox also 
contains a Dealunay algorithm that uses the centroids to calculate 
contiguous entities.

I'm no GIS expert (I'm an applied econometrician) and the code I've written 
seems to work. The Distance package also works with my real data which 
are the centroids of the counties in Connecticut and I tested it with 
Euclidean, Cityblock, and SqEuclidean.

Once again, thank you all for the help. I'm hoping this thread can help 
others and get them into Julia as well.

Regards,
Don

On Thursday, July 3, 2014 12:34:28 AM UTC-4, Patrick O'Leary wrote:

 (/me removes Julia hat, replaces with navigation hat)

 That distance measurement probably isn't going to do quite what you want 
 it to do, particularly over such a large area. It's an okay approximation 
 in CONUS, but biased--1 degree latitude isn't the same distance as 1 degree 
 longitude. For points further north this gets worse. You have a few options 
 for more accuracy:

 You can put all of your coordinates in a local level frame. All the 
 necessary formulas and constants are conveniently provided in MATLAB's 
 Aerospace Blockset documentation: 
 http://www.mathworks.com/help/aeroblks/llatoflatearth.html. Choose 
 something like the centroid lat/lon as the initial values, which will act 
 as the origin of the reference frame. Note that lat/lon will need to be 
 converted to radians.

 For longer distances, you will get better results with something like the 
 haversine formula, which accounts for the curvature of the earth (
 https://en.wikipedia.org/wiki/Haversine_formula). That would probably 
 make for a decent contribution to Distance.jl, and is pretty easy to 
 implement.

 For (substantial!) extra credit, you could implement a version that 
 accounts for the flattening of the ellipsoid (
 https://en.wikipedia.org/wiki/Geodesics_on_an_ellipsoid, 
 https://en.wikipedia.org/wiki/Vincenty%27s_formulae). That would probably 
 not be accepted by Distance.jl, and is unlikely to be worthwhile for your 
 application.

 On Wednesday, July 2, 2014 7:05:21 PM UTC-5, Donald Lacombe wrote:

 Actually, another issue with the type information.

 Instead of using random coordinates, I used actual latitude and longitude 
 coordinates from a CSV file like so:

 data = readcsv(ct_coord.csv);

 temp = [xc';yc']

 2x8 Array{Any,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


 Now, the type here is {Any,2} and this seems to be giving the Distance 
 package problems:



 R = pairwise(Euclidean(),temp)

 MethodError(fptype,(Any,))


 Julia seems to be unhappy with this type but I'm not sure how to correct 
 this which would seem to require that the data being read in are of Float 
 type.


 Again, thank you for all of the help.


 Don



 On Wednesday, July 2, 2014 7:10:09 PM UTC-4, Donald Lacombe wrote:

 Patrick,

 Thank you for the reply! I will try your suggestion and hopefully get a 
 better grasp of these issues.

 Thank you all gain for all of the help. It cannot be stressed enough 
 that the friendly and truly helpful comments will help the Julia community 
 grow.

 Don

 On Wednesday, July 2, 2014 5:31:55 PM UTC-4, Patrick O'Leary wrote:

 The second parameter to the Array type is the number of dimensions of 
 the array. Row vectors are implemented with 2-dimensional arrays 
 (a.k.a., 
 matrices), since Vectors (a.k.a. Array{T, 1} for all types T) are 
 interpreted as columns.

 Since pairwise() generates all of the distances simultaneously, you'll 
 need to slice its results to use it in sortperm(), following your original 
 loop implementation. Each slice will need to be a Vector. If you slice 
 columns, you get Vector slices for free. If you slice rows, the shape will 
 be maintained, so you get a row back. You can flatten this to a Vector 
 with 
 the vec() function.

 Hope this helps!

 On Wednesday, July 2, 2014 4:23:52 PM UTC-5, Donald Lacombe wrote:

 After thinking about this issue, another couple of question came to 
 mind.

 When I create a row vector, say xc = randn(1,5), the results indicate 
 that the vector is e.g.

 xc = rand(1,5)

 1x5 Array{Float64,2}:

  0.989421  0.718148  0.499914  0.13024  0.939856


 The type is {Float64,2}.


 Questions:


 1) What does the 2 mean in {Float64,2}? If I create xc = rand(5) the 
 type is {Float64,1}.


 2) The sortperm function only seems to accept {Float64,1} types. The 
 Euclidean distance code create types of {Float64,2} but I need 
 {Float64,1} for the sortperm function. Do I need to convert these before 
 using this command?


 I'm probably missing something really basic so apologies for the query.


 Thanks,

 Don


 On Wednesday, July 2, 2014 4:20:55 PM UTC-4, Donald Lacombe wrote:

 Greetings, part 2...

 After doing some tests, I found

[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-03 Thread Donald Lacombe
Patrick (and others),

Another issue that has reared it's ugly head is that when I read the data 
using the Data Frames package, I get the following:

data = readtable(ct_coord_2.csv,header=false)

8x2 DataFrame:

  x1  x2

[1,]-73.3712  41.225

[2,]-72.1065 41.4667

[3,]-73.2453 41.7925

[4,]-71.9876   41.83

[5,]-72.3365  41.855

[6,]-72.7328 41.8064

[7,]-72.5231 41.4354

[8,]-72.8999 41.3488


julia xc = data[:,1]

8-element DataArray{Float64,1}:

 -73.3712

 -72.1065

 -73.2453

 -71.9876

 -72.3365

 -72.7328

 -72.5231

 -72.8999


julia yc = data[:,2]

8-element DataArray{Float64,1}:

 41.225 

 41.4667

 41.7925

 41.83  

 41.855 

 41.8064

 41.4354

 41.3488


julia xc=xc'

1x8 DataArray{Float64,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999


julia yc=yc'

1x8 DataArray{Float64,2}:

 41.225  41.4667  41.7925  41.83  41.855  41.8064  41.4354  41.3488


julia temp = [xc;yc]

2x8 DataArray{Float64,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

  41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


julia R = pairwise(Euclidean(),temp)

MethodError(At_mul_B!,(

8x8 Array{Float64,2}:

 2.7273e-316   2.7273e-316   2.67478e-315  …  2.7273e-316   2.7273e-316 

 2.67736e-315  2.67736e-315  2.67736e-315 2.72726e-316  2.72726e-316

 2.67727e-315  2.67727e-315  2.67727e-315 2.67727e-315  2.67727e-315

 2.67727e-315  2.67727e-315  2.67727e-315 2.67727e-315  2.67727e-315

 4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324

 2.76235e-318  2.76235e-318  2.76235e-318  …  2.76235e-318  2.76235e-318

 4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324

 4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324,


2x8 DataArray{Float64,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

  41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488,


2x8 DataArray{Float64,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

  41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488))


I do not think that the Distance package likes the types that is input into the 
function, i.e. the vectors are DataArrays instead of Arrays. It works just fine 
when I used Tony's idea:


 julia data = readcsv(ct_coord_2.csv,Float64)

8x2 Array{Float64,2}:

 -73.3712  41.225 

 -72.1065  41.4667

 -73.2453  41.7925

 -71.9876  41.83  

 -72.3365  41.855 

 -72.7328  41.8064

 -72.5231  41.4354

 -72.8999  41.3488


julia xc = data[:,1]

8-element Array{Float64,1}:

 -73.3712

 -72.1065

 -73.2453

 -71.9876

 -72.3365

 -72.7328

 -72.5231

 -72.8999


julia yc = data[:,2]

8-element Array{Float64,1}:

 41.225 

 41.4667

 41.7925

 41.83  

 41.855 

 41.8064

 41.4354

 41.3488


julia xc=xc'

1x8 Array{Float64,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999


julia yc=yc'

1x8 Array{Float64,2}:

 41.225  41.4667  41.7925  41.83  41.855  41.8064  41.4354  41.3488


julia temp = [xc;yc]

2x8 Array{Float64,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

  41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


julia R = pairwise(Euclidean(),temp)

8x8 Array{Float64,2}:

 0.0   1.28762   0.581327  1.51014   …  0.863479  0.873799  0.487347

 1.28762   0.0   1.18451   0.382214 0.712542  0.417808  0.802085

 0.581327  1.18451   0.0   1.25833  0.512668  0.805673  0.562309

 1.51014   0.382214  1.25833   0.0  0.745667  0.665227  1.03141 

 1.21144   0.451294  0.910982  0.349837 0.399323  0.459258  0.757372

 0.863479  0.712542  0.512668  0.745667  …  0.0   0.426208  0.487124

 0.873799  0.417808  0.805673  0.665227 0.426208  0.0   0.386557

 0.487347  0.802085  0.562309  1.03141  0.487124  0.386557  0.0 


There seems to be some issue with the Distance package not accepting Data 
Frames. Of course, the readcsv works fine but this might be an issue for others 
as well.


Thanks,

Don



On Thursday, July 3, 2014 6:49:18 PM UTC-4, Patrick O'Leary wrote:

 On Thursday, July 3, 2014 5:36:23 PM UTC-5, Donald Lacombe wrote:

 I'm no GIS expert (I'm an applied econometrician) and the code I've 
 written seems to work. The Distance package also works with my real data 
 which are the centroids of the counties in Connecticut and I tested it with 
 Euclidean, Cityblock, and SqEuclidean.


 Glad you got something working. Whether those distances are accurate 
 enough depends on how the points are arranged and what you plan to do with 
 it--I can see where it wouldn't make much difference in this case. I can't 
 let the statisticians and image processing folks have all the technical 
 conversation fun in this mailing list, though!



[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-03 Thread Donald Lacombe
Dear Johan,

Tank you for this bit of information. It seems to me that the DataFrames 
package might become the defacto manner in which data is imported into 
Julia but I'd still like to have both options open.

Thank you all again,
Don 

On Thursday, July 3, 2014 7:54:49 PM UTC-4, Johan Sigfrids wrote:

 You can use dropna() to convert a DataArray to a Array. This will 
 obviously drop any missing values. 

 On Friday, July 4, 2014 2:08:55 AM UTC+3, Donald Lacombe wrote:

 Patrick (and others),

 Another issue that has reared it's ugly head is that when I read the data 
 using the Data Frames package, I get the following:

 data = readtable(ct_coord_2.csv,header=false)

 8x2 DataFrame:

   x1  x2

 [1,]-73.3712  41.225

 [2,]-72.1065 41.4667

 [3,]-73.2453 41.7925

 [4,]-71.9876   41.83

 [5,]-72.3365  41.855

 [6,]-72.7328 41.8064

 [7,]-72.5231 41.4354

 [8,]-72.8999 41.3488


 julia xc = data[:,1]

 8-element DataArray{Float64,1}:

  -73.3712

  -72.1065

  -73.2453

  -71.9876

  -72.3365

  -72.7328

  -72.5231

  -72.8999


 julia yc = data[:,2]

 8-element DataArray{Float64,1}:

  41.225 

  41.4667

  41.7925

  41.83  

  41.855 

  41.8064

  41.4354

  41.3488


 julia xc=xc'

 1x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999


 julia yc=yc'

 1x8 DataArray{Float64,2}:

  41.225  41.4667  41.7925  41.83  41.855  41.8064  41.4354  41.3488


 julia temp = [xc;yc]

 2x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


 julia R = pairwise(Euclidean(),temp)

 MethodError(At_mul_B!,(

 8x8 Array{Float64,2}:

  2.7273e-316   2.7273e-316   2.67478e-315  …  2.7273e-316   2.7273e-316 

  2.67736e-315  2.67736e-315  2.67736e-315 2.72726e-316  2.72726e-316

  2.67727e-315  2.67727e-315  2.67727e-315 2.67727e-315  2.67727e-315

  2.67727e-315  2.67727e-315  2.67727e-315 2.67727e-315  2.67727e-315

  4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324

  2.76235e-318  2.76235e-318  2.76235e-318  …  2.76235e-318  2.76235e-318

  4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324

  4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324,


 2x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488,


 2x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488))


 I do not think that the Distance package likes the types that is input into 
 the function, i.e. the vectors are DataArrays instead of Arrays. It works 
 just fine when I used Tony's idea:


  julia data = readcsv(ct_coord_2.csv,Float64)

 8x2 Array{Float64,2}:

  -73.3712  41.225 

  -72.1065  41.4667

  -73.2453  41.7925

  -71.9876  41.83  

  -72.3365  41.855 

  -72.7328  41.8064

  -72.5231  41.4354

  -72.8999  41.3488


 julia xc = data[:,1]

 8-element Array{Float64,1}:

  -73.3712

  -72.1065

  -73.2453

  -71.9876

  -72.3365

  -72.7328

  -72.5231

  -72.8999


 julia yc = data[:,2]

 8-element Array{Float64,1}:

  41.225 

  41.4667

  41.7925

  41.83  

  41.855 

  41.8064

  41.4354

  41.3488


 julia xc=xc'

 1x8 Array{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999


 julia yc=yc'

 1x8 Array{Float64,2}:

  41.225  41.4667  41.7925  41.83  41.855  41.8064  41.4354  41.3488


 julia temp = [xc;yc]

 2x8 Array{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


 julia R = pairwise(Euclidean(),temp)

 8x8 Array{Float64,2}:

  0.0   1.28762   0.581327  1.51014   …  0.863479  0.873799  0.487347

  1.28762   0.0   1.18451   0.382214 0.712542  0.417808  0.802085

  0.581327  1.18451   0.0   1.25833  0.512668  0.805673  0.562309

  1.51014   0.382214  1.25833   0.0  0.745667  0.665227  1.03141 

  1.21144   0.451294  0.910982  0.349837 0.399323  0.459258  0.757372

  0.863479  0.712542  0.512668  0.745667  …  0.0   0.426208  0.487124

  0.873799  0.417808  0.805673  0.665227 0.426208  0.0   0.386557

  0.487347  0.802085  0.562309  1.03141  0.487124  0.386557  0.0 


 There seems to be some issue with the Distance package not accepting Data 
 Frames. Of course, the readcsv works fine but this might be an issue for 
 others as well.


 Thanks,

 Don



 On Thursday, July 3, 2014 6:49:18 PM UTC-4, Patrick O'Leary wrote:

 On Thursday, July 3, 2014 5:36:23 PM UTC-5, Donald Lacombe wrote:

 I'm no GIS expert (I'm an applied econometrician) and the code I've 
 written seems to work. The Distance package also works with my

[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-03 Thread Donald Lacombe
Johan,

I think there may be an issue with the Data Frames package as I get the 
following:

julia data = readtable(test.csv,header=false)

6x2 DataFrame:

x1 x2

[1,] 1  7

[2,] 2  8

[3,] 3  9

[4,] 4 10

[5,] 5 11

[6,] 6 12



julia convert(Array,data)

MethodError(convert,(Array{T,N},6x2 DataFrame:

x1 x2

[1,] 1  7

[2,] 2  8

[3,] 3  9

[4,] 4 10

[5,] 5 11

[6,] 6 12

))


julia dropna(data)

ErrorException(dropna not defined)


I read the documentation and they both say the same thing but it doesn't seem 
to work in my case.


Thoughts?


Thanks,

Don


On Thursday, July 3, 2014 7:54:49 PM UTC-4, Johan Sigfrids wrote:

 You can use dropna() to convert a DataArray to a Array. This will 
 obviously drop any missing values. 

 On Friday, July 4, 2014 2:08:55 AM UTC+3, Donald Lacombe wrote:

 Patrick (and others),

 Another issue that has reared it's ugly head is that when I read the data 
 using the Data Frames package, I get the following:

 data = readtable(ct_coord_2.csv,header=false)

 8x2 DataFrame:

   x1  x2

 [1,]-73.3712  41.225

 [2,]-72.1065 41.4667

 [3,]-73.2453 41.7925

 [4,]-71.9876   41.83

 [5,]-72.3365  41.855

 [6,]-72.7328 41.8064

 [7,]-72.5231 41.4354

 [8,]-72.8999 41.3488


 julia xc = data[:,1]

 8-element DataArray{Float64,1}:

  -73.3712

  -72.1065

  -73.2453

  -71.9876

  -72.3365

  -72.7328

  -72.5231

  -72.8999


 julia yc = data[:,2]

 8-element DataArray{Float64,1}:

  41.225 

  41.4667

  41.7925

  41.83  

  41.855 

  41.8064

  41.4354

  41.3488


 julia xc=xc'

 1x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999


 julia yc=yc'

 1x8 DataArray{Float64,2}:

  41.225  41.4667  41.7925  41.83  41.855  41.8064  41.4354  41.3488


 julia temp = [xc;yc]

 2x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


 julia R = pairwise(Euclidean(),temp)

 MethodError(At_mul_B!,(

 8x8 Array{Float64,2}:

  2.7273e-316   2.7273e-316   2.67478e-315  …  2.7273e-316   2.7273e-316 

  2.67736e-315  2.67736e-315  2.67736e-315 2.72726e-316  2.72726e-316

  2.67727e-315  2.67727e-315  2.67727e-315 2.67727e-315  2.67727e-315

  2.67727e-315  2.67727e-315  2.67727e-315 2.67727e-315  2.67727e-315

  4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324

  2.76235e-318  2.76235e-318  2.76235e-318  …  2.76235e-318  2.76235e-318

  4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324

  4.94066e-324  4.94066e-324  4.94066e-324 9.88131e-324  4.94066e-324,


 2x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488,


 2x8 DataArray{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488))


 I do not think that the Distance package likes the types that is input into 
 the function, i.e. the vectors are DataArrays instead of Arrays. It works 
 just fine when I used Tony's idea:


  julia data = readcsv(ct_coord_2.csv,Float64)

 8x2 Array{Float64,2}:

  -73.3712  41.225 

  -72.1065  41.4667

  -73.2453  41.7925

  -71.9876  41.83  

  -72.3365  41.855 

  -72.7328  41.8064

  -72.5231  41.4354

  -72.8999  41.3488


 julia xc = data[:,1]

 8-element Array{Float64,1}:

  -73.3712

  -72.1065

  -73.2453

  -71.9876

  -72.3365

  -72.7328

  -72.5231

  -72.8999


 julia yc = data[:,2]

 8-element Array{Float64,1}:

  41.225 

  41.4667

  41.7925

  41.83  

  41.855 

  41.8064

  41.4354

  41.3488


 julia xc=xc'

 1x8 Array{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999


 julia yc=yc'

 1x8 Array{Float64,2}:

  41.225  41.4667  41.7925  41.83  41.855  41.8064  41.4354  41.3488


 julia temp = [xc;yc]

 2x8 Array{Float64,2}:

  -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

   41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


 julia R = pairwise(Euclidean(),temp)

 8x8 Array{Float64,2}:

  0.0   1.28762   0.581327  1.51014   …  0.863479  0.873799  0.487347

  1.28762   0.0   1.18451   0.382214 0.712542  0.417808  0.802085

  0.581327  1.18451   0.0   1.25833  0.512668  0.805673  0.562309

  1.51014   0.382214  1.25833   0.0  0.745667  0.665227  1.03141 

  1.21144   0.451294  0.910982  0.349837 0.399323  0.459258  0.757372

  0.863479  0.712542  0.512668  0.745667  …  0.0   0.426208  0.487124

  0.873799  0.417808  0.805673  0.665227 0.426208  0.0   0.386557

  0.487347  0.802085  0.562309  1.03141  0.487124  0.386557  0.0 


 There seems to be some issue with the Distance package

[julia-users] Using Distance Package for knn Weight Matrix

2014-07-02 Thread Donald Lacombe
Greetings!

I am a new Julia user and have the following issue. I am writing code to 
calculate a knn spatial weight matrix to estimate a spatial econometric 
model. Using the MATLAB code from Jim LeSage's Econometrics Toolbox, I 
converted that code into Julia as follows:

xc = rand(8);

yc = rand(8);

n = length(xc);

k = 2;

distance = zeros(n);

nnlist = zeros(n,k);

tempw = zeros(n,n);


for i=1:n;

xi = xc[i,1];

yi = yc[i,1];

distance = (xc - xi*ones(n)).^2 + (yc - yi*ones(n)).^2

temp = sortperm(distance)

nnlist[i,1:k] = temp[2:k+1,1]';

end


for i=1:n

tempw[i,nnlist[i,:]] = 1;

end


W = tempw/k;


This is a toy example and I was wondering if I can use the Distance package 
to simplify the distance = (xc - xi*ones(n)).^2 + (yc - yi*ones(n)).^2 formula. 
I tried using the pairwise option like so:


R = pairwise(Euclidean(),xc,yc) but received the following message:


R = pairwise(Euclidean(),xc,yc)

MethodError(pairwise,(Euclidean(),[0.05961066617957589,0.018538084399339905,0.39282193332224646,0.7006919213133509,0.5099836895629475,0.8448415935222402,0.2985674570217043,0.8022287058003177],[0.5808687231553928,0.9655167324458858,0.026306556019434435,0.6565373244339141,0.11927452074471412,0.11873635450496622,0.6271632933770979,0.7081439899673692]))

I'd like to be able to utilize the Distance package but am a bit stumped. The 
code as written works but it's bugging me that I cannot seem to get the above 
command to work. I also get the following error when loading Distance:


using Distance

Warning: could not import Base.foldl into NumericExtensions

Warning: could not import Base.foldr into NumericExtensions

Warning: could not import Base.sum! into NumericExtensions

Warning: could not import Base.maximum! into NumericExtensions

Warning: could not import Base.minimum! into NumericExtensions


If there is anyone who can address this I'd greatly appreciate it. 
Incidentally, the help on this group is one reason I am making the change.


Don






[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-02 Thread Donald Lacombe
Greetings, part 2...

After doing some tests, I found that this code works just fine:

julia using Distance

Warning: could not import Base.foldl into NumericExtensions

Warning: could not import Base.foldr into NumericExtensions

Warning: could not import Base.sum! into NumericExtensions

Warning: could not import Base.maximum! into NumericExtensions

Warning: could not import Base.minimum! into NumericExtensions



julia xc = rand(8)'

1x8 Array{Float64,2}:

 0.30107  0.479169  0.230607  0.65126  0.386921  0.455316  0.921496  0.244873


julia yc = rand(8)'

1x8 Array{Float64,2}:

 0.866199  0.767794  0.76163  0.262925  …  0.742765  0.980952  0.424966


julia temp = [xc ; yc]

2x8 Array{Float64,2}:

 0.30107   0.479169  0.230607  0.65126   …  0.455316  0.921496  0.244873

 0.866199  0.767794  0.76163   0.262925 0.742765  0.980952  0.424966


julia R = pairwise(Euclidean(),temp)

8x8 Array{Float64,2}:

 0.0   0.203477   0.126094  0.697548  …  0.197554   0.630949  0.444797

 0.203477  0.00.248638  0.533393 0.0345751  0.491009  0.415241

 0.126094  0.248638   0.0   0.652423 0.225499   0.724865  0.336966

 0.697548  0.533393   0.652423  0.0  0.518306   0.767196  0.437502

 0.376201  0.283308   0.304834  0.355027 0.252287   0.719136  0.160613

 0.197554  0.0345751  0.225499  0.518306  …  0.00.523505  0.381159

 0.630949  0.491009   0.724865  0.767196 0.523505   0.0   0.87575 

 0.444797  0.415241   0.336966  0.437502 0.381159   0.87575   0.0 


As expected, the diagonal elements are 0 and everything seems to work fine 
although I did not test other distance formulas.

I'm posing this so others can see the results and possibly make use of this 
themselves.

Is there a central repository for user submitted solutions or examples? If 
so I'd like to post this so others can see or alternatively someone else is 
free to use this (and the original code I posted) for future reference.

Thank you all again for your help! I'm going to be coding some models this 
summer and will probably be back to bug everyone again.

Thanks,
Don

On Wednesday, July 2, 2014 10:22:58 AM UTC-4, Donald Lacombe wrote:

 Greetings!

 I am a new Julia user and have the following issue. I am writing code to 
 calculate a knn spatial weight matrix to estimate a spatial econometric 
 model. Using the MATLAB code from Jim LeSage's Econometrics Toolbox, I 
 converted that code into Julia as follows:

 xc = rand(8);

 yc = rand(8);

 n = length(xc);

 k = 2;

 distance = zeros(n);

 nnlist = zeros(n,k);

 tempw = zeros(n,n);


 for i=1:n;

 xi = xc[i,1];

 yi = yc[i,1];

 distance = (xc - xi*ones(n)).^2 + (yc - yi*ones(n)).^2

 temp = sortperm(distance)

 nnlist[i,1:k] = temp[2:k+1,1]';

 end


 for i=1:n

 tempw[i,nnlist[i,:]] = 1;

 end


 W = tempw/k;


 This is a toy example and I was wondering if I can use the Distance package 
 to simplify the distance = (xc - xi*ones(n)).^2 + (yc - yi*ones(n)).^2 
 formula. I tried using the pairwise option like so:


 R = pairwise(Euclidean(),xc,yc) but received the following message:


 R = pairwise(Euclidean(),xc,yc)

 MethodError(pairwise,(Euclidean(),[0.05961066617957589,0.018538084399339905,0.39282193332224646,0.7006919213133509,0.5099836895629475,0.8448415935222402,0.2985674570217043,0.8022287058003177],[0.5808687231553928,0.9655167324458858,0.026306556019434435,0.6565373244339141,0.11927452074471412,0.11873635450496622,0.6271632933770979,0.7081439899673692]))

 I'd like to be able to utilize the Distance package but am a bit stumped. The 
 code as written works but it's bugging me that I cannot seem to get the above 
 command to work. I also get the following error when loading Distance:


 using Distance

 Warning: could not import Base.foldl into NumericExtensions

 Warning: could not import Base.foldr into NumericExtensions

 Warning: could not import Base.sum! into NumericExtensions

 Warning: could not import Base.maximum! into NumericExtensions

 Warning: could not import Base.minimum! into NumericExtensions


 If there is anyone who can address this I'd greatly appreciate it. 
 Incidentally, the help on this group is one reason I am making the change.


 Don






[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-02 Thread Donald Lacombe
After thinking about this issue, another couple of question came to mind.

When I create a row vector, say xc = randn(1,5), the results indicate that 
the vector is e.g.

xc = rand(1,5)

1x5 Array{Float64,2}:

 0.989421  0.718148  0.499914  0.13024  0.939856


The type is {Float64,2}.


Questions:


1) What does the 2 mean in {Float64,2}? If I create xc = rand(5) the type is 
{Float64,1}.


2) The sortperm function only seems to accept {Float64,1} types. The Euclidean 
distance code create types of {Float64,2} but I need {Float64,1} for the 
sortperm function. Do I need to convert these before using this command?


I'm probably missing something really basic so apologies for the query.


Thanks,

Don


On Wednesday, July 2, 2014 4:20:55 PM UTC-4, Donald Lacombe wrote:

 Greetings, part 2...

 After doing some tests, I found that this code works just fine:

 julia using Distance

 Warning: could not import Base.foldl into NumericExtensions

 Warning: could not import Base.foldr into NumericExtensions

 Warning: could not import Base.sum! into NumericExtensions

 Warning: could not import Base.maximum! into NumericExtensions

 Warning: could not import Base.minimum! into NumericExtensions



 julia xc = rand(8)'

 1x8 Array{Float64,2}:

  0.30107  0.479169  0.230607  0.65126  0.386921  0.455316  0.921496  0.244873


 julia yc = rand(8)'

 1x8 Array{Float64,2}:

  0.866199  0.767794  0.76163  0.262925  …  0.742765  0.980952  0.424966


 julia temp = [xc ; yc]

 2x8 Array{Float64,2}:

  0.30107   0.479169  0.230607  0.65126   …  0.455316  0.921496  0.244873

  0.866199  0.767794  0.76163   0.262925 0.742765  0.980952  0.424966


 julia R = pairwise(Euclidean(),temp)

 8x8 Array{Float64,2}:

  0.0   0.203477   0.126094  0.697548  …  0.197554   0.630949  0.444797

  0.203477  0.00.248638  0.533393 0.0345751  0.491009  0.415241

  0.126094  0.248638   0.0   0.652423 0.225499   0.724865  0.336966

  0.697548  0.533393   0.652423  0.0  0.518306   0.767196  0.437502

  0.376201  0.283308   0.304834  0.355027 0.252287   0.719136  0.160613

  0.197554  0.0345751  0.225499  0.518306  …  0.00.523505  0.381159

  0.630949  0.491009   0.724865  0.767196 0.523505   0.0   0.87575 

  0.444797  0.415241   0.336966  0.437502 0.381159   0.87575   0.0 


 As expected, the diagonal elements are 0 and everything seems to work fine 
 although I did not test other distance formulas.

 I'm posing this so others can see the results and possibly make use of 
 this themselves.

 Is there a central repository for user submitted solutions or examples? If 
 so I'd like to post this so others can see or alternatively someone else is 
 free to use this (and the original code I posted) for future reference.

 Thank you all again for your help! I'm going to be coding some models this 
 summer and will probably be back to bug everyone again.

 Thanks,
 Don

 On Wednesday, July 2, 2014 10:22:58 AM UTC-4, Donald Lacombe wrote:

 Greetings!

 I am a new Julia user and have the following issue. I am writing code to 
 calculate a knn spatial weight matrix to estimate a spatial econometric 
 model. Using the MATLAB code from Jim LeSage's Econometrics Toolbox, I 
 converted that code into Julia as follows:

 xc = rand(8);

 yc = rand(8);

 n = length(xc);

 k = 2;

 distance = zeros(n);

 nnlist = zeros(n,k);

 tempw = zeros(n,n);


 for i=1:n;

 xi = xc[i,1];

 yi = yc[i,1];

 distance = (xc - xi*ones(n)).^2 + (yc - yi*ones(n)).^2

 temp = sortperm(distance)

 nnlist[i,1:k] = temp[2:k+1,1]';

 end


 for i=1:n

 tempw[i,nnlist[i,:]] = 1;

 end


 W = tempw/k;


 This is a toy example and I was wondering if I can use the Distance 
 package to simplify the distance = (xc - xi*ones(n)).^2 + (yc - 
 yi*ones(n)).^2 formula. I tried using the pairwise option like so:


 R = pairwise(Euclidean(),xc,yc) but received the following message:


 R = pairwise(Euclidean(),xc,yc)

 MethodError(pairwise,(Euclidean(),[0.05961066617957589,0.018538084399339905,0.39282193332224646,0.7006919213133509,0.5099836895629475,0.8448415935222402,0.2985674570217043,0.8022287058003177],[0.5808687231553928,0.9655167324458858,0.026306556019434435,0.6565373244339141,0.11927452074471412,0.11873635450496622,0.6271632933770979,0.7081439899673692]))

 I'd like to be able to utilize the Distance package but am a bit stumped. 
 The code as written works but it's bugging me that I cannot seem to get the 
 above command to work. I also get the following error when loading Distance:


 using Distance

 Warning: could not import Base.foldl into NumericExtensions

 Warning: could not import Base.foldr into NumericExtensions

 Warning: could not import Base.sum! into NumericExtensions

 Warning: could not import Base.maximum! into NumericExtensions

 Warning: could not import Base.minimum! into NumericExtensions


 If there is anyone who can address this I'd greatly appreciate it. 
 Incidentally, the help

[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-02 Thread Donald Lacombe
Patrick,

Thank you for the reply! I will try your suggestion and hopefully get a 
better grasp of these issues.

Thank you all gain for all of the help. It cannot be stressed enough that 
the friendly and truly helpful comments will help the Julia community grow.

Don

On Wednesday, July 2, 2014 5:31:55 PM UTC-4, Patrick O'Leary wrote:

 The second parameter to the Array type is the number of dimensions of the 
 array. Row vectors are implemented with 2-dimensional arrays (a.k.a., 
 matrices), since Vectors (a.k.a. Array{T, 1} for all types T) are 
 interpreted as columns.

 Since pairwise() generates all of the distances simultaneously, you'll 
 need to slice its results to use it in sortperm(), following your original 
 loop implementation. Each slice will need to be a Vector. If you slice 
 columns, you get Vector slices for free. If you slice rows, the shape will 
 be maintained, so you get a row back. You can flatten this to a Vector with 
 the vec() function.

 Hope this helps!

 On Wednesday, July 2, 2014 4:23:52 PM UTC-5, Donald Lacombe wrote:

 After thinking about this issue, another couple of question came to mind.

 When I create a row vector, say xc = randn(1,5), the results indicate 
 that the vector is e.g.

 xc = rand(1,5)

 1x5 Array{Float64,2}:

  0.989421  0.718148  0.499914  0.13024  0.939856


 The type is {Float64,2}.


 Questions:


 1) What does the 2 mean in {Float64,2}? If I create xc = rand(5) the type 
 is {Float64,1}.


 2) The sortperm function only seems to accept {Float64,1} types. The 
 Euclidean distance code create types of {Float64,2} but I need {Float64,1} 
 for the sortperm function. Do I need to convert these before using this 
 command?


 I'm probably missing something really basic so apologies for the query.


 Thanks,

 Don


 On Wednesday, July 2, 2014 4:20:55 PM UTC-4, Donald Lacombe wrote:

 Greetings, part 2...

 After doing some tests, I found that this code works just fine:

 julia using Distance

 Warning: could not import Base.foldl into NumericExtensions

 Warning: could not import Base.foldr into NumericExtensions

 Warning: could not import Base.sum! into NumericExtensions

 Warning: could not import Base.maximum! into NumericExtensions

 Warning: could not import Base.minimum! into NumericExtensions



 julia xc = rand(8)'

 1x8 Array{Float64,2}:

  0.30107  0.479169  0.230607  0.65126  0.386921  0.455316  0.921496  
 0.244873


 julia yc = rand(8)'

 1x8 Array{Float64,2}:

  0.866199  0.767794  0.76163  0.262925  …  0.742765  0.980952  0.424966


 julia temp = [xc ; yc]

 2x8 Array{Float64,2}:

  0.30107   0.479169  0.230607  0.65126   …  0.455316  0.921496  0.244873

  0.866199  0.767794  0.76163   0.262925 0.742765  0.980952  0.424966


 julia R = pairwise(Euclidean(),temp)

 8x8 Array{Float64,2}:

  0.0   0.203477   0.126094  0.697548  …  0.197554   0.630949  0.444797

  0.203477  0.00.248638  0.533393 0.0345751  0.491009  0.415241

  0.126094  0.248638   0.0   0.652423 0.225499   0.724865  0.336966

  0.697548  0.533393   0.652423  0.0  0.518306   0.767196  0.437502

  0.376201  0.283308   0.304834  0.355027 0.252287   0.719136  0.160613

  0.197554  0.0345751  0.225499  0.518306  …  0.00.523505  0.381159

  0.630949  0.491009   0.724865  0.767196 0.523505   0.0   0.87575 

  0.444797  0.415241   0.336966  0.437502 0.381159   0.87575   0.0 


 As expected, the diagonal elements are 0 and everything seems to work 
 fine although I did not test other distance formulas.

 I'm posing this so others can see the results and possibly make use of 
 this themselves.

 Is there a central repository for user submitted solutions or examples? 
 If so I'd like to post this so others can see or alternatively someone else 
 is free to use this (and the original code I posted) for future reference.

 Thank you all again for your help! I'm going to be coding some models 
 this summer and will probably be back to bug everyone again.

 Thanks,
 Don

 On Wednesday, July 2, 2014 10:22:58 AM UTC-4, Donald Lacombe wrote:

 Greetings!

 I am a new Julia user and have the following issue. I am writing code 
 to calculate a knn spatial weight matrix to estimate a spatial econometric 
 model. Using the MATLAB code from Jim LeSage's Econometrics Toolbox, I 
 converted that code into Julia as follows:

 xc = rand(8);

 yc = rand(8);

 n = length(xc);

 k = 2;

 distance = zeros(n);

 nnlist = zeros(n,k);

 tempw = zeros(n,n);


 for i=1:n;

 xi = xc[i,1];

 yi = yc[i,1];

 distance = (xc - xi*ones(n)).^2 + (yc - yi*ones(n)).^2

 temp = sortperm(distance)

 nnlist[i,1:k] = temp[2:k+1,1]';

 end


 for i=1:n

 tempw[i,nnlist[i,:]] = 1;

 end


 W = tempw/k;


 This is a toy example and I was wondering if I can use the Distance 
 package to simplify the distance = (xc - xi*ones(n)).^2 + (yc - 
 yi*ones(n)).^2 formula. I tried using the pairwise option like so:


 R = pairwise(Euclidean(),xc,yc

[julia-users] Re: Using Distance Package for knn Weight Matrix

2014-07-02 Thread Donald Lacombe
Actually, another issue with the type information.

Instead of using random coordinates, I used actual latitude and longitude 
coordinates from a CSV file like so:

data = readcsv(ct_coord.csv);

temp = [xc';yc']

2x8 Array{Any,2}:

 -73.3712  -72.1065  -73.2453  -71.9876  …  -72.7328  -72.5231  -72.8999

  41.22541.4667   41.7925   41.8341.8064   41.4354   41.3488


Now, the type here is {Any,2} and this seems to be giving the Distance package 
problems:



R = pairwise(Euclidean(),temp)

MethodError(fptype,(Any,))


Julia seems to be unhappy with this type but I'm not sure how to correct this 
which would seem to require that the data being read in are of Float type.


Again, thank you for all of the help.


Don



On Wednesday, July 2, 2014 7:10:09 PM UTC-4, Donald Lacombe wrote:

 Patrick,

 Thank you for the reply! I will try your suggestion and hopefully get a 
 better grasp of these issues.

 Thank you all gain for all of the help. It cannot be stressed enough that 
 the friendly and truly helpful comments will help the Julia community grow.

 Don

 On Wednesday, July 2, 2014 5:31:55 PM UTC-4, Patrick O'Leary wrote:

 The second parameter to the Array type is the number of dimensions of the 
 array. Row vectors are implemented with 2-dimensional arrays (a.k.a., 
 matrices), since Vectors (a.k.a. Array{T, 1} for all types T) are 
 interpreted as columns.

 Since pairwise() generates all of the distances simultaneously, you'll 
 need to slice its results to use it in sortperm(), following your original 
 loop implementation. Each slice will need to be a Vector. If you slice 
 columns, you get Vector slices for free. If you slice rows, the shape will 
 be maintained, so you get a row back. You can flatten this to a Vector with 
 the vec() function.

 Hope this helps!

 On Wednesday, July 2, 2014 4:23:52 PM UTC-5, Donald Lacombe wrote:

 After thinking about this issue, another couple of question came to mind.

 When I create a row vector, say xc = randn(1,5), the results indicate 
 that the vector is e.g.

 xc = rand(1,5)

 1x5 Array{Float64,2}:

  0.989421  0.718148  0.499914  0.13024  0.939856


 The type is {Float64,2}.


 Questions:


 1) What does the 2 mean in {Float64,2}? If I create xc = rand(5) the type 
 is {Float64,1}.


 2) The sortperm function only seems to accept {Float64,1} types. The 
 Euclidean distance code create types of {Float64,2} but I need {Float64,1} 
 for the sortperm function. Do I need to convert these before using this 
 command?


 I'm probably missing something really basic so apologies for the query.


 Thanks,

 Don


 On Wednesday, July 2, 2014 4:20:55 PM UTC-4, Donald Lacombe wrote:

 Greetings, part 2...

 After doing some tests, I found that this code works just fine:

 julia using Distance

 Warning: could not import Base.foldl into NumericExtensions

 Warning: could not import Base.foldr into NumericExtensions

 Warning: could not import Base.sum! into NumericExtensions

 Warning: could not import Base.maximum! into NumericExtensions

 Warning: could not import Base.minimum! into NumericExtensions



 julia xc = rand(8)'

 1x8 Array{Float64,2}:

  0.30107  0.479169  0.230607  0.65126  0.386921  0.455316  0.921496  
 0.244873


 julia yc = rand(8)'

 1x8 Array{Float64,2}:

  0.866199  0.767794  0.76163  0.262925  …  0.742765  0.980952  0.424966


 julia temp = [xc ; yc]

 2x8 Array{Float64,2}:

  0.30107   0.479169  0.230607  0.65126   …  0.455316  0.921496  0.244873

  0.866199  0.767794  0.76163   0.262925 0.742765  0.980952  0.424966


 julia R = pairwise(Euclidean(),temp)

 8x8 Array{Float64,2}:

  0.0   0.203477   0.126094  0.697548  …  0.197554   0.630949  0.444797

  0.203477  0.00.248638  0.533393 0.0345751  0.491009  0.415241

  0.126094  0.248638   0.0   0.652423 0.225499   0.724865  0.336966

  0.697548  0.533393   0.652423  0.0  0.518306   0.767196  0.437502

  0.376201  0.283308   0.304834  0.355027 0.252287   0.719136  0.160613

  0.197554  0.0345751  0.225499  0.518306  …  0.00.523505  0.381159

  0.630949  0.491009   0.724865  0.767196 0.523505   0.0   0.87575 

  0.444797  0.415241   0.336966  0.437502 0.381159   0.87575   0.0 


 As expected, the diagonal elements are 0 and everything seems to work 
 fine although I did not test other distance formulas.

 I'm posing this so others can see the results and possibly make use of 
 this themselves.

 Is there a central repository for user submitted solutions or examples? 
 If so I'd like to post this so others can see or alternatively someone 
 else 
 is free to use this (and the original code I posted) for future reference.

 Thank you all again for your help! I'm going to be coding some models 
 this summer and will probably be back to bug everyone again.

 Thanks,
 Don

 On Wednesday, July 2, 2014 10:22:58 AM UTC-4, Donald Lacombe wrote:

 Greetings!

 I am a new Julia user and have the following issue. I am writing code