Mailman 3 May 2018 - NumPy-Discussion

mixed mode arithmetic
by Neal Becker July 11, 2023

July 11, 2023

I've been browsing the numpy source. I'm wondering about mixed-mode arithmetic on arrays. I believe the way numpy handles this is that it never does mixed arithmetic, but instead converts arrays to a common type. Arguably, that might be efficient for a mix of say, double and float. Maybe not. But for a mix of complex and a scalar type (say, CDouble * Double), it's clearly suboptimal in efficiency. So, do I understand this correctly? If so, is that something we should improve?

4 6

Invalid value encoutered : how to prevent numpy.where to do this?
by Eric Emsellem Feb. 18, 2023

Feb. 18, 2023

Dear all, I have a code using lots of "numpy.where" to make some constrained calculations as in: data = arange(10) result = np.where(data == 0, 0., 1./data) # or data1 = arange(10) data2 = arange(10)+1.0 result = np.where(data1 > data2, np.sqrt(data1-data2), np.sqrt(data2-data2)) which then produces warnings like: /usr/bin/ipython:1: RuntimeWarning: invalid value encountered in sqrt or for the first example: /usr/bin/ipython:1: RuntimeWarning: divide by zero encountered in divide How … [View More]

5 4

Re: [Numpy-discussion] example reading binary Fortran file
by Neil Martinsen-Burrell July 22, 2021

July 22, 2021

David Froger <david.froger.info <at> gmail.com> writes: > Hy,My question is about reading Fortran binary file (oh no this question > again...) I've posted this before, but I finally got it cleaned up for the Cookbook. For this purpose I use a subclass of file that has methods for reading unformatted Fortran data. See http://www.scipy.org/Cookbook/FortranIO/FortranFile. I'd gladly see this in numpy or scipy somewhere, but I'm not sure where it belongs. > program … [View More]makeArray > implicit none > integer,parameter:: nx=10,ny=20 > real(4),dimension(nx,ny):: ux,uy,p > integer :: i,j > open(11,file='uxuyp.bin',form='unformatted') > do i = 1,nx > do j = 1,ny > ux(i,j) = real(i*j) > uy(i,j) = real(i)/real(j) > p (i,j) = real(i) + real(j) > enddo > enddo > write(11) ux,uy > write(11) p > close(11) > end program makeArray When I run the above program compiled with gfortran on my Intel Mac, I can read it back with:: >>> import numpy as np >>> from fortranfile import FortranFile >>> f=FortranFile('uxuyp.bin', endian='<') >>> uxuy = f.readReals(prec='f') # 'f' for default reals >>> len(uxuy) 400 >>> ux = np.array(uxuy[:200]).reshape((20,10)).T >>> uy = np.array(uxuy[200:]).reshape((20,10)).T >>> p = f.readReals('f').reshape((20,10)).T >>> ux array([[ 1., 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20.], [ 2., 4., 6., 8., 10., 12., 14., 16., 18., 20., 22., 24., 26., 28., 30., 32., 34., 36., 38., 40.], [ 3., 6., 9., 12., 15., 18., 21., 24., 27., 30., 33., 36., 39., 42., 45., 48., 51., 54., 57., 60.], [ 4., 8., 12., 16., 20., 24., 28., 32., 36., 40., 44., 48., 52., 56., 60., 64., 68., 72., 76., 80.], [ 5., 10., 15., 20., 25., 30., 35., 40., 45., 50., 55., 60., 65., 70., 75., 80., 85., 90., 95., 100.], [ 6., 12., 18., 24., 30., 36., 42., 48., 54., 60., 66., 72., 78., 84., 90., 96., 102., 108., 114., 120.], [ 7.,Proxy-Connection: keep-alive Cache-Control: max-age=0 14., 21., 28., 35., 42., 49., 56., 63., 70., 77., 84., 91., 98., 105., 112., 119., 126., 133., 140.], [ 8., 16., 24., 32., 40., 48., 56., 64., 72., 80., 88., 96., 104., 112., 120., 128., 136., 144., 152., 160.], [ 9., 18., 27., 36., 45., 54., 63., 72., 81., 90., 99., 108., 117., 126., 135., 144., 153., 162., 171., 180.], [ 10., 20., 30., 40., 50., 60., 70., 80., 90., 100., 110., 120., 130., 140., 150., 160., 170., 180., 190., 200.]]) >>> uy array([[ 1. , 0.5 , 0.33333334, 0.25 , 0.2 , 0.16666667, 0.14285715, 0.125 , 0.11111111, 0.1 , 0.09090909, 0.08333334, 0.07692308, 0.07142857, 0.06666667, 0.0625 , 0.05882353, 0.05555556, 0.05263158, 0.05 ], [ 2. , 1. , 0.66666669, 0.5 , 0.40000001, 0.33333334, 0.2857143 , 0.25 , 0.22222222, 0.2 , 0.18181819, 0.16666667, 0.15384616, 0.14285715, 0.13333334, 0.125 , 0.11764706, 0.11111111, 0.10526316, 0.1 ], [ 3. , 1.5 , 1. , 0.75 , 0.60000002, 0.5 , 0.42857143, 0.375 , 0.33333334, 0.30000001, 0.27272728, 0.25 , 0.23076923, 0.21428572, 0.2 , 0.1875 , 0.17647059, 0.16666667, 0.15789473, 0.15000001], [ 4. , 2. , 1.33333337, 1. , 0.80000001, 0.66666669, 0.5714286 , 0.5 , 0.44444445, 0.40000001, 0.36363637, 0.33333334, 0.30769232, 0.2857143 , 0.26666668, 0.25 , 0.23529412, 0.22222222, 0.21052632, 0.2 ], [ 5. , 2.5 , 1.66666663, 1.25 , 1. , 0.83333331, 0.71428573, 0.625 , 0.55555558, 0.5 , 0.45454547, 0.41666666, 0.38461539, 0.35714287, 0.33333334, 0.3125 , 0.29411766, 0.27777779, 0.2631579 , 0.25 ], [ 6. , 3. , 2. , 1.5 , 1.20000005, 1. , 0.85714287, 0.75 , 0.66666669, 0.60000002, 0.54545456, 0.5 , 0.46153846, 0.42857143, 0.40000001, 0.375 , 0.35294119, 0.33333334, 0.31578946, 0.30000001], [ 7. , 3.5 , 2.33333325, 1.75 , 1.39999998, 1.16666663, 1. , 0.875 , 0.77777779, 0.69999999, 0.63636363, 0.58333331, 0.53846157, 0.5 , 0.46666667, 0.4375 , 0.41176471, 0.3888889 , 0.36842105, 0.34999999], [ 8. , 4. , 2.66666675, 2. , 1.60000002, 1.33333337, 1.14285719, 1. , 0.8888889 , 0.80000001, 0.72727275, 0.66666669, 0.61538464, 0.5714286 , 0.53333336, 0.5 , 0.47058824, 0.44444445, 0.42105263, 0.40000001], [ 9. , 4.5 , 3. , 2.25 , 1.79999995, 1.5 , 1.28571427, 1.125 , 1. , 0.89999998, 0.81818181, 0.75 , 0.69230771, 0.64285713, 0.60000002, 0.5625 , 0.52941179, 0.5 , 0.47368422, 0.44999999], [ 10. , 5. , 3.33333325, 2.5 , 2. , 1.66666663, 1.42857146, 1.25 , 1.11111116, 1. , 0.90909094, 0.83333331, 0.76923078, 0.71428573, 0.66666669, 0.625 , 0.58823532, 0.55555558, 0.52631581, 0.5 ]]) >>> p array([[ 2., 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21.], [ 3., 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22.], [ 4., 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23.], [ 5., 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24.], [ 6., 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25.], [ 7., 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26.], [ 8., 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27.], [ 9., 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28.], [ 10., 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28., 29.], [ 11., 12., 13., 14., 15., 16., 17., 18., 19., 20., 21., 22., 23., 24., 25., 26., 27., 28., 29., 30.]]) Note that you have to provide the shape information for ux and uy because fortran writes them together as a stream of 400 numbers. -Neil [View Less]

10 10

Allowing broadcasting of code dimensions in generalized ufuncs
by Marten van Kerkwijk July 3, 2018

July 3, 2018

Hi All, Following on a PR combining the ability to provide fixed and flexible dimensions [1] (useful for, e.g., 3-vector input with a signature like `(3),(3)->(3)`, and for `matmul`, resp.; based on earlier PRs by Jaime [2] and Matt (Picus) [3]), I've now made a PR with a further enhancement, which allows one can indicate that a core dimension can be broadcast [4]. A particular use case is `all_equal`, a new function suggested in a stalled PR by Matt (Harrigan) [5], which compares two … [View More]

10 42

A roadmap for NumPy - longer term planning
by Matti Picus June 4, 2018

June 4, 2018

At the recent NumPy sprint at BIDS (thanks to those who made the trip) we spent some time brainstorming about a roadmap for NumPy, in the spirit of similar work that was done for Jupyter. The idea is that a document with wide community acceptance can guide the work of the full-time developer(s), and be a source of ideas for expanding development efforts. I put the document up at https://github.com/numpy/numpy/wiki/NumPy-Roadmap, and hope to discuss it at a BOF session during SciPy in … [View More]

12 14

Where to discuss NEPs (was: Re: new NEP: np.AbstractArray and np.asabstractarray)
by Nathaniel Smith May 31, 2018

May 31, 2018

On Thu, Mar 8, 2018 at 7:06 AM, Marten van Kerkwijk <m.h.vankerkwijk(a)gmail.com> wrote: > Hi Nathaniel, > > Overall, hugely in favour! For detailed comments, it would be good to > have a link to a PR; could you put that up? Well, there's a PR here: https://github.com/numpy/numpy/pull/10706 But, this raises a question :-). (One which also came up here: https://github.com/numpy/numpy/pull/10704#issuecomment-371684170) There are sensible two workflows we could use (or at … [View More]

7 15

Inconsistent results for the covariance matrix between scipy.optimize.curve_fit and numpy.polyfit
by Jonathan Tammo Siebert May 29, 2018

May 29, 2018

Hi, I hope this is the appropriate place to ask something like this, otherwise please let me know (or feel free to ignore this). Also I hope that I do not misunderstood something or did some silly mistake. If so, please let me know as well! TLDR: When scaling the covariance matrix based on the residuals, scipy.optimize.curve_fit uses a factor of chisq(popt)/(M-N) (with M=number of point, N=number of parameters) and numpy.polyfit uses chisq(popt)/(M-N-2). I am wondering which is correct. I am … [View More]somewhat confused about different results I am getting for the covariance matrix of a simple linear fit, when comparing `scipy.optimize.curve_fit` and `numpy.polyfit`. I am aware, that `curve_fit` solves the more general non-linear problem numerically, while `polyfit` finds an analytical solution to the linear problem. However, both converge to the same solution, so I suspect that this difference is not important here. The difference, I am curious about is not in the returned parameters but in the estimate of the corresponding covariance matrix. As I understand, there are two different ways to estimate it, based either on the absolute values of the provided uncertainties or by interpreting those only as weights and then scaling the matrix to produce an appropriate reduced chisq. To that end, curve_fit has the parameter: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.cur ve_fit.html: "absolute_sigma : bool, optional If True, sigma is used in an absolute sense and the estimated parameter covariance pcov reflects these absolute values. If False, only the relative magnitudes of the sigma values matter. The returned parameter covariance matrix pcov is based on scaling sigma by a constant factor. This constant is set by demanding that the reduced chisq for the optimal parameters popt when using the scaled sigma equals unity. In other words, sigma is scaled to match the sample variance of the residuals after the fit. Mathematically, pcov(absolute_sigma=False) = pcov(absolute_sigma=True) * chisq(popt)/(M-N)" https://docs.scipy.org/doc/numpy/reference/generated/numpy.polyfit.html on the other hand, does not say anything about how the covariance matrix is estimated. To my understanding, its default should correspond to `absolute_sigma=False` for `curve_fit`. As `polyfit` has a weight parameter instead of an uncertainty parameter, I guess the difference in default behavior is not that surprising. However, even when specifying `absolute_sigma=False`, `curve_fit` and `polyfit` produce different covariance matrices as the applied scaling factors are chisq(popt)/(M-N-2) for `polyfit` (https://github.com/numpy/numpy/blob/6a58e25703cbecb6786faa09a04ae2ec82 21348b/numpy/lib/polynomial.py#L598-L605) and chisq(popt)/(M-N) for `curve_fit` (https://github.com/scipy/scipy/blob/607a21e07dad234f8e63fcf03b7994137a 3ccd5b/scipy/optimize/minpack.py#L781-L782). The argument given in a comment to the scaling `polyfit` is: "Some literature ignores the extra -2.0 factor in the denominator, but it is included here because the covariance of Multivariate Student-T (which is implied by a Bayesian uncertainty analysis) includes it. Plus, it gives a slightly more conservative estimate of uncertainty.", but honestly, in a quick search, I was not able to find any literature not ignoring the extra "factor". But obviously, I could very well be misunderstanding something. Nonetheless, as `curve_fit` ignores it as well, I was wondering whether those two shouldn't give consistent results and if so, which would be the correct solution. Best, Jonathan [View Less]

2 3

matmul as a ufunc
by Matti Picus May 29, 2018

May 29, 2018

I have made progress with resolving the issue that matmul, the operation which implements `a @ b`, is not a ufunc [2]. Discussion on the issue, which prevents the __array_ufunc__ mechanism for overriding matmul on subclasses of ndarray, yeilded two approaches: - create a wrapper that can convince the ufunc mechanism to call __array_ufunc__ even on functions that are not true ufuncs - expand the semantics of core signatures so that a single matmul ufunc can implement matrix-matrix, vector-… [View More]

7 9

ANN: EuroScipy 2018
by Valerio Maggio May 29, 2018

May 29, 2018

*** Apologies if you receive multiple copies *** Dear Colleagues, We are delighted to invite you to join us for the **11th European Conference on Python in Science**. The EuroSciPy 2018 (https://www.euroscipy.org/2018/) Conference will be organised by Fondazione Bruno Kessler (FBK) and will take place from August 28 to September 1 in **Trento, Italy**. The EuroSciPy meeting is a cross-disciplinary gathering focused on the use and development of the Python language in scientific research. … [View More]This event strives to bring together both users and developers of scientific tools, as well as academic research and state of the art industry. The conference will be structured as it follows: Aug, 28-29 : Tutorials and Hands-on Aug, 30-31 : Main Conference Sep, 1 : Sprint TOPICS OF INTEREST: =================== Presentations of scientific tools and libraries using the Python language, including but not limited to: - Algorithms implemented or exposed in Python - Astronomy - Data Visualisation - Deep Learning & AI - Earth, Ocean and Geo Science - General-purpose Python tools that can be of special interest to the scientific community. - Image Processing - Materials Science - Parallel computing - Political and Social Sciences - Project Jupyter - Reports on the use of Python in scientific achievements or ongoing projects. - Robotics & IoT - Scientific data flow and persistence - Scientific visualization - Simulation - Statistics - Vector and array manipulation - Web applications and portals for science and engineering - 3D Printing CALL FOR PROPOSALS: =================== EuroScipy will accept three different kinds of contributions: * Regular Talks: standard talks for oral presentations, allocated in time slots of 15, or 30 minutes, depending on your preference and scheduling constraints. Each time slot considers a Q&A session at the end of the talk (at least, 5 mins). * Hands-on Tutorials: These are beginner or advanced training sessions to dive into the subject with all details. These sessions are 90 minutes long, and the audience will be strongly encouraged to bring a laptop to experiment. * Poster: EuroScipy will host two poster sessions during the two days of Main Conference. So attendees and students are highly encourage to present their work and/or preliminary results as posters. Proposals should be submitted using the EuroScipy submission system at https://pretalx.com/euroscipy18. **Submission deadline is May, 31st 2018**. REGISTRATION & FEES: ==================== To register to EuroScipy 2018, please go to http://euroscipy2018.eventbrite.co.uk or to http://www.euroscipy.org/2018/ Fees: ----- | Tutorials | Student* | Academic/Individual | Industry| |-----------------------------|----------|---------------------|----------| | Early Bird (till July, 1st) | € 50,00 | € 70,00 | € 125,00 | | Regular (till Aug, 5th | € 100,00 | € 110,00 | € 250,00 | | Late (till Aug, 22nd) | € 135,00 | € 135,00 | € 300,00 | | Main Conference | Student* | Academic/Individual | Industry| |-----------------------------|----------|---------------------|----------| | Early Bird (till July, 1st) | € 50,00 | € 70,00 | € 125,00 | | Regular (till Aug, 5th | € 100,00 | € 110,00 | € 250,00 | | Late (till Aug, 22nd) | € 135,00 | € 135,00 | € 300,00 | * A proof of student status will be required at time of the registration. kind regards, Valerio EuroScipy 2018 Organising Committee, Email: info(a)euroscipy.org Website: http://www.euroscipy.org/2018 twitter: @euroscipy [View Less]

1 0

Adding take_along_axis and put_along_axis functions
by Eric Wieser May 28, 2018

May 28, 2018

These functions provide a vectorized way of using one array to look up items in another. In particular, they extend the 1d: a = np.array([4, 5, 6, 1, 2, 3]) b = np.array(["four", "five", "six", "one", "two", "three"]) i = a.argsort() b_sorted = b[i] To work for higher-dimensions: a = np.array([[4, 1], [5, 2], [6, 3]]) b = np.array([["four", "one"], ["five", "two"], ["six", "three"]]) i = a.argsort(axis=1) b_sorted = np.take_along_axis(b, i, axis=1) put_along_axis is the obvious but less … [View More]

2 1